apps130311_water.mapcs_plot

import sys
import pickle as pickle

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap

import mapinmedata
import irfpy.util.gridsphere

def main(filename):
    if filename.endswith('_matrix.txt'):
        print('Reading data %s' % filename)
        fp = open(filename)
        nobs = pickle.load(fp)
        sumdata = pickle.load(fp)
        md = irfpy.util.gridsphere.SimpleGridSphere()
        lonarr, latarr = md.get_cgrid()
    else:
        print('Reading data %s' % filename)
        m = mapinmedata.orbitdata_reader(filename)
        print('...done')
    
        md = m  # For alias
        lonarr, latarr = md.get_cgrid()
        nobs = md.get_counts()
        print('Max ndata=', nobs.max())

        datalist = md.datalist
        sundata = [np.array(d).sum() for d in datalist]
        sumdata = np.array(sundata)
        sumdata.shape=(180, 360)
        fp2 = open(filename+'_matrix.txt', 'w')
        pickle.dump(nobs, fp2)
        pickle.dump(sumdata, fp2)
        fp2.close()

    avg = sumdata / nobs
    fig = plt.figure(figsize=(12, 6))
    ax1 = fig.add_subplot(1, 2, 1)
    ax2 = fig.add_subplot(1, 2, 2)

    m = Basemap(projection='spaeqd', boundinglat=-40, lon_0=90, ax=ax1)
#    m.drawcoastlines(color='0.7')
    m.drawparallels(np.arange(-80, 81, 10))
    m.drawmeridians(np.arange(-180, 181, 30), labels=[1, 0, 0, 1])
    x, y = m(lonarr, latarr)
    nobs = np.ma.masked_less(nobs, 0.1)
    if not nobs.mask.all():
                m.pcolor(x, y, nobs)
    ax1.set_title('Number of observation')
        
    m = Basemap(projection='spaeqd', boundinglat=-40, lon_0=90, ax=ax2)
#    m.drawcoastlines(color='0.7')
    m.drawparallels(np.arange(-80, 81, 10))
    m.drawmeridians(np.arange(-180, 181, 30), labels=[0, 1, 0, 1])
    x, y = m(lonarr, latarr)
    if not nobs.mask.all():
        m.pcolor(x, y, avg, vmax=np.percentile(avg, 99.7))
    ax2.set_title('Simple sum of CENA counts')

    fig.savefig(filename+'.png')

if __name__ == '__main__':
    main(sys.argv[1])