apps120712_spice1205.mission_summary_monthly

''' Mission summary from 1205 spice kernel.

.. image:: ../../../src/scripts/apps120712_spice1205/mission_summary_203001.png
    :width: 100%
'''

import datetime

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

import pyana.pep.juice_spice1205 as jspice1205
import pyana.juice.jspice as jspice
jspice.init()

def main():

    missiontable = jspice1205.JuiceSummary.table

    # Mission JOI.
    tjoi = missiontable['G1'][1]
    print('JOI is', tjoi)
    t0 = datetime.datetime(tjoi.year, tjoi.month, 1)

    # Mission end
    tend = missiontable['END'][1]
    print('END is', tend)
    t1 = tend


    ### Full orbit
    # Create time list  (from tjoi to tend with 1 day tick)
    resolution = datetime.timedelta(hours=12)
    tlist_full = [mpldates.num2date(t).replace(tzinfo=None) for t in mpldates.drange(t0, t1, resolution)]

    ### Position in JSE frame
    poslist_full = []
    for t in tlist_full:
        poslist_full.append(jspice.get_position(t))
    poslist_full = np.array(poslist_full)

    cmap = matplotlib.cm.get_cmap('jet', 6)
    vmin = 0
    vmax = 18

    ### Monthly plots

    for yr in (2030, 2031, 2032, 2033):
        for mo in range(1, 13):
            tm0 = datetime.datetime(yr, mo, 1)    
            if tm0 > t1: continue
            if tm0 < t0: tm0 = t0

            if mo != 12:
                tm1 = datetime.datetime(yr, mo + 1, 1)
            else:
                tm1 = datetime.datetime(yr + 1, 1, 1)
            if tm1 < t0: continue
            if tm1 > t1: tm1 = t1

            resolution = datetime.timedelta(hours=1)
            tlist_monthly = [mpldates.num2date(t).replace(tzinfo=None)
                            for t in mpldates.drange(tm0, tm1, resolution)]
            ### Position in JSE frame
            poslist_monthly = []
            for t in tlist_monthly:
                poslist_monthly.append(jspice.get_position(t))
            poslist_monthly = np.array(poslist_monthly)
            print(np.abs(poslist_monthly[:, 2]).max() / 71500)

            ### Velocity
            vellist_monthly = []
            for t in tlist_monthly:
                vellist_monthly.append(jspice.get_velocity(t))
            vellist_monthly = np.array(vellist_monthly)
            vellist_monthly = np.sqrt((vellist_monthly ** 2).sum(1))
#            print vellist_monthly.max()

            ### Moon pos
            moon1pos = []
            moon2pos = []
            moon3pos = []
            moon4pos = []
            for t in tlist_monthly[::3]:
                moon1pos.append(jspice.get_position(t, target='IO'))
                moon2pos.append(jspice.get_position(t, target='EUROPA'))
                moon3pos.append(jspice.get_position(t, target='GANYMEDE'))
                moon4pos.append(jspice.get_position(t, target='CALLISTO'))
            moon1pos = np.array(moon1pos)
            moon2pos = np.array(moon2pos)
            moon3pos = np.array(moon3pos)
            moon4pos = np.array(moon4pos)
            

            rj = 71500

            plt.figure(figsize=(16, 8))

            plt.subplot(131)

            plt.plot(poslist_full[:, 0] / rj, poslist_full[:, 1] / rj, ':', color='gray')
            plt.plot(moon1pos[:, 0] / rj,
                     moon1pos[:, 1] / rj, 'b-')
            plt.plot(moon2pos[:, 0] / rj,
                     moon2pos[:, 1] / rj, 'b-')
            plt.plot(moon3pos[:, 0] / rj,
                     moon3pos[:, 1] / rj, 'b-')
            plt.plot(moon4pos[:, 0] / rj,
                     moon4pos[:, 1] / rj, 'b-')
                    
#            plt.plot(poslist_monthly[:, 0] / rj, poslist_monthly[:, 1] / rj, 'r-', lw=3)
            plt.scatter(poslist_monthly[:, 0] / rj, poslist_monthly[:, 1] / rj, c=vellist_monthly, edgecolor='none', vmin=vmin, vmax=vmax, cmap=cmap)
            plt.gca().set_aspect('equal')
            plt.title('S/C position [Rj]')
            plt.xlabel('X (JSE) [Rj]')
            plt.ylabel('Y (JSE) [Rj]')

            plt.subplot(132)

            plt.plot(poslist_full[:, 0] / rj, poslist_full[:, 1] / rj, ':', color='gray')
            plt.plot(moon1pos[:, 0] / rj,
                     moon1pos[:, 1] / rj, 'b-')
            plt.plot(moon2pos[:, 0] / rj,
                     moon2pos[:, 1] / rj, 'b-')
            plt.plot(moon3pos[:, 0] / rj,
                     moon3pos[:, 1] / rj, 'b-')
            plt.plot(moon4pos[:, 0] / rj,
                     moon4pos[:, 1] / rj, 'b-')
                    
#            plt.plot(poslist_monthly[:, 0] / rj, poslist_monthly[:, 1] / rj, 'r-', lw=3)
            sc = plt.scatter(poslist_monthly[:, 0] / rj, poslist_monthly[:, 1] / rj, c=vellist_monthly, edgecolor='none', vmin=vmin, vmax=vmax, cmap=cmap)
            plt.xlim(-40, 30)
            plt.ylim(-40, 30)
            plt.gca().set_aspect('equal')
            plt.title('S/C position [Rj]')
            plt.xlabel('X (JSE) [Rj]')
            plt.ylabel('Y (JSE) [Rj]')

            plt.subplot(133)

            plt.plot(poslist_full[:, 0] / rj, poslist_full[:, 2] / rj, ':', color='gray')
            plt.plot(moon1pos[:, 0] / rj,
                     moon1pos[:, 2] / rj, 'b-')
            plt.plot(moon2pos[:, 0] / rj,
                     moon2pos[:, 2] / rj, 'b-')
            plt.plot(moon3pos[:, 0] / rj,
                     moon3pos[:, 2] / rj, 'b-')
            plt.plot(moon4pos[:, 0] / rj,
                     moon4pos[:, 2] / rj, 'b-')

            sc = plt.scatter(poslist_monthly[:, 0] / rj, poslist_monthly[:, 2] / rj, c=vellist_monthly, edgecolor='none', vmin=vmin, vmax=vmax, cmap=cmap)
            plt.xlim(-40, 30)
            plt.ylim(-40, 30)
            plt.gca().set_aspect('equal')
            plt.title('S/C position [Rj]')
            plt.xlabel('X (JSE) [Rj]')
            plt.ylabel('Z (JSE) [Rj]')
            


            cb = plt.colorbar(sc)
            cb.set_label('Velocity of S/C [km/s]')

            plt.savefig('mission_summary_%02d%02d.png' % (yr, mo))

    return


if __name__ == "__main__":
    main()