tessellatedsphere_searcharea

''' A sample to plot searching tessellated data
'''
import time

import numpy as np
from matplotlib.mlab import griddata
import matplotlib.pyplot as plt
import matplotlib.cm
from matplotlib.patches import Polygon
from matplotlib.collections import PatchCollection


from irfpy.util import tessellatedsphere

def main():
    '''Main script'''

    ### Here you instance tessellated sphere of 4 dimension.
    ts = tessellatedsphere.RecurrentTessellatedUnitSphere(4)
    ts.createLinkList()
    
    def pole_easthem_func(lons, lats):
        return np.logical_and(lons >= 0, np.abs(lats) > 66.7)


    fig = plt.figure()

    for istr, strategy in enumerate(['array_full', 'single_full', 'single_link', 'array_link']):

        print('Stragety: {}'.format(strategy))
        _t0 = time.time()
        if strategy.endswith('full'):
            neidx = ts.search_area(pole_easthem_func, strategy=strategy)
        else:
            neidx = ts.search_area(pole_easthem_func, [15.0, 83.0], strategy=strategy)
        print('    ... ellapsed = {} s'.format(time.time() - _t0))

        pcoll = []
        ccoll = []   # Saving color.

        for idx, tri in enumerate(ts.tess_tris):
            vert = np.array(tri.getVertices())
            lonv = np.arctan2(vert[:, 1], vert[:, 0]) * 180. / np.pi
            latv = np.arcsin(vert[:, 2]) * 180. / np.pi

            lonlat = np.array([lonv, latv])    # (2, 3) array

            if lonv.max() - lonv.min() >= 180:  # Too wide, then ignore.
                continue

            if idx in neidx:
                isfill = True
                ccoll.append(1)
            else:
                ccoll.append(0)
            polygon = Polygon(lonlat.T, closed=True)
            pcoll.append(polygon)

        ax = fig.add_subplot(2, 2, istr + 1)
        p = PatchCollection(pcoll, cmap=matplotlib.cm.jet, edgecolor='none')
        p.set_array(np.array(ccoll))

        ax.add_collection(p)

        ax.set_aspect('equal')

        ax.set_xlim(-180, 180)
        ax.set_ylim(-90, 90)
        ax.set_title(strategy)
    
    fig.savefig('tessellatedsphere_searcharea.png')

if __name__ == "__main__":
    main()