apps160222_jdc.fovplotΒΆ

To look at the FOV under the S/C lon-lat frame

""" To look at the FOV under the S/C lon-lat frame
"""
import numpy as np
import matplotlib.pyplot as plt

def frame_plot_bins(ax=None, tilt=0):
    """

    :param ax: Axis to be plotted
    :param tilt: Tilt angle (in degrees)
    :return: axis
    """
    if ax is None:
        fig = plt.figure()
        ax = fig.add_subplot(111)

    ### FOV of JDC
    from irfpy.jdc import fov1
    from irfpy.jdc import frame1

    for ielev in range(4):
        for iazim in range(16):

            ### Convert from index to angle
            elev_jdc = fov1.elev_pix_center(ielev)
            azim_jdc = fov1.azim_pix_center(iazim)
            while azim_jdc > 180:
                azim_jdc -= 360
            while azim_jdc <= -180:
                azim_jdc += 360

            ### Convert angle to vector (in JDC frame)
            jdc_vec = frame1.angles2jdc(elev_jdc, azim_jdc)

            ### Convert JDC vector to SC frame
            sc_vec = frame1.jdc2nsc(jdc_vec, tilt=tilt)


            elev_sc = np.rad2deg(np.arcsin(sc_vec[2]))
            azim_sc = np.rad2deg(np.arctan2(sc_vec[1], sc_vec[0]))

            print(ielev, iazim, elev_jdc, azim_jdc, jdc_vec, sc_vec, elev_sc, azim_sc)

            ax.plot([azim_sc], [elev_sc], 'bo')

    frame_plot_blockage(ax)

    ### Decoration
    ax.set_xlabel('Longitude S/C')
    ax.set_ylabel('Latitude S/C')

    ax.set_xlim(-180, 180)
    ax.set_ylim(-90, 90)

    return ax

def frame_plot_blockage(ax=None):
    if ax is None:
        fig = plt.figure()
        ax = fig.add_subplot(111)

    ### Satellite blockage
    import matplotlib.patches as mplp

    sc = mplp.Rectangle((-180, 0), 90, 90, fc='y')
    ax.add_patch(sc)

    sc = mplp.Rectangle((90, 0), 90, 90, fc='y')
    ax.add_patch(sc)

    nozzle = mplp.Rectangle((-180, -30), 10, 30, fc='y')
    ax.add_patch(nozzle)

    nozzle = mplp.Rectangle((170, -30), 10, 30, fc='y')
    ax.add_patch(nozzle)

    return ax

def frame_plot_limit(ax=None, tilt=0, **argv):
    """

    :param ax: Axis to be plotted
    :param tilt: Tilt angle (in degrees)
    :return: axis
    """
    if ax is None:
        fig = plt.figure()
        ax = fig.add_subplot(111)

    ### FOV of JDC
    from irfpy.jdc import fov1
    from irfpy.jdc import frame1

    ### lower and upper boundary.
    for elev_jdc in [90, ]:
        for azim_jdc in range(-180, 181):
            ### Convert angle to vector (in JDC frame)
            jdc_vec = frame1.angles2jdc(elev_jdc, azim_jdc)

            ### Convert JDC vector to SC frame
            sc_vec = frame1.jdc2nsc(jdc_vec, tilt=tilt)


            elev_sc = np.rad2deg(np.arcsin(sc_vec[2]))
            azim_sc = np.rad2deg(np.arctan2(sc_vec[1], sc_vec[0]))

#            print(ielev, iazim, elev_jdc, azim_jdc, jdc_vec, sc_vec, elev_sc, azim_sc)

            ax.plot([azim_sc], [elev_sc], 'r.', **argv)


    frame_plot_blockage(ax)

    ### Decoration
    ax.set_xlabel('Longitude S/C')
    ax.set_ylabel('Latitude S/C')

    ax.set_xlim(-180, 180)
    ax.set_ylim(-90, 90)

    return ax


def main():
    """ Plot will be a map in S/C frame.
    """
    ### Something should be plotted here
    fig = plt.figure(figsize=(20, 10))
    ax = fig.add_subplot(221)
    ax = frame_plot_bins(ax=ax, tilt=0)
    ax = frame_plot_limit(ax=ax, tilt=0)
    ax.set_title('Tilt=0 deg')

    ax = fig.add_subplot(222)
    ax = frame_plot_bins(ax=ax, tilt=15)
    ax = frame_plot_limit(ax=ax, tilt=15)
    ax.set_title('Tilt=15 deg')

    ax = fig.add_subplot(223)
    ax = frame_plot_bins(ax=ax, tilt=30)
    ax = frame_plot_limit(ax=ax, tilt=30)
    ax.set_title('Tilt=30 deg')

    ax = fig.add_subplot(224)
    ax = frame_plot_bins(ax=ax, tilt=45)
    ax = frame_plot_limit(ax=ax, tilt=45)
    ax.set_title('Tilt=45 deg')

    ### Plot
    fig.savefig('fovplot.png')


if __name__ == '__main__':
    main()