apps121023_jesframe.sundir

''' To check Sun direction
'''

import datetime
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.dates as mpldates

import pyana.juice.jspice as js
js.init()

def main():
    '''Main script'''
    
    fig = plt.figure(figsize=(12, 6))
    ax = fig.add_subplot(111)

    t0 = datetime.datetime(2011, 1, 1)
    t0m = mpldates.date2num(t0)
    t1 = datetime.datetime(2031, 1, 31)
    dt = datetime.timedelta(minutes=1440)

    t = t0

    xax = np.array([1, 0, 0])
    jselist = []
    jeqslist = []
    tlist = []

    while t <= t1:
        
        sundir_jse = js.get_position(t, frame='JSE', target='Sun', origin='Jupiter')
        sundir_jse = np.array(sundir_jse)
        sundir_jse = sundir_jse / np.sqrt((sundir_jse * sundir_jse).sum())   # Unit vector
        sundir_angle = np.rad2deg(np.arccos((xax * sundir_jse).sum()))    # Angle from X-axis
        jselist.append(sundir_angle)

        sundir_jeqs = js.get_position(t, frame='JEQS', target='Sun', origin='Jupiter')
        sundir_jeqs = np.array(sundir_jeqs)
        sundir_jeqs = sundir_jeqs / np.sqrt((sundir_jeqs * sundir_jeqs).sum())   # Unit vector
        sundir_angle = np.rad2deg(np.arccos((xax * sundir_jeqs).sum()))    # Angle from X-axis
        jeqslist.append(sundir_angle)

        tlist.append(mpldates.date2num(t))

        t = t + dt

    ax.plot_date(tlist, jselist, 'b', label='JSE X-Sun angle')
    ax.plot_date(tlist, jeqslist, 'r', label='JEqS X-Sun angle')
    ax.set_yscale('log')
    ax.legend(loc='best')
    fig.autofmt_xdate()

    fig.savefig('sundir.png')

if __name__ == "__main__":
    main()