#!/usr/bin/env python

from pyhdf.SD import SD, SDC 

import numpy as np
import matplotlib.pyplot as plt
import matplotlib as mpl
import matplotlib.cm as cm
import my_sivo

#----------------------------------------------------------------------------------------#
#
# Inputs
#
#----------------------------------------------------------------------------------------#

myd06_file_name = './inputs/broken_clouds/MYD06_L2.A2007219.2010.006.2014053202546.hdf'

#----------------------------------------------------------------------------------------#
#
# Read HDF Files
#
#----------------------------------------------------------------------------------------#

myd06 = SD(myd06_file_name, SDC.READ)

myd06_cloud_phase_optical_properties_c5 = myd06.select('cloud_top_temperature_1km')

myd06_cloud_phase_optical_properties_c5_data = myd06_cloud_phase_optical_properties_c5.get()

print "Read HDF done"

#----------------------------------------------------------------------------------------#
#
# Read HDF Files
#
#----------------------------------------------------------------------------------------#


myd06_cloud_phase_optical_properties_c5_data = (myd06_cloud_phase_optical_properties_c5_data + 15000) * 0.01

print myd06_cloud_phase_optical_properties_c5_data

'''
print myd06_cloud_phase_optical_properties_c5_data.shape
data_shape = myd06_cloud_phase_optical_properties_c5_data.shape

for i in range(data_shape[0]):
	for j in range(data_shape[1]):
		if myd06_cloud_phase_optical_properties_c5_data[i,j] > 285.0:
			print myd06_cloud_phase_optical_properties_c5_data[i,j], i,j
'''

#cmap = [(0.75,0.75,0.75)] + [(my_sivo.nasa (i)) for i in xrange(1,256)] 
cmap = [(0.0,0.0,0.0)] + [(my_sivo.nasa (i)) for i in xrange(1,256)] 
cmap = mpl.colors.ListedColormap(cmap)
bounds = [100,200,205,210,215,220,225,230,235,240,245,250,255,260,265,270,275,280,285,290,300]
norm = mpl.colors.BoundaryNorm(bounds, cmap.N)

#img = plt.imshow(np.fliplr(myd06_cloud_phase_optical_properties_c5_data), cmap=cmap, norm=norm, interpolation='nearest', origin='lower', vmin =200.0, vmax=300.0)
img = plt.imshow(np.fliplr(myd06_cloud_phase_optical_properties_c5_data), cmap=cmap, norm=norm, interpolation='none', origin='lower')


bounds = [200,205,210,215,220,225,230,235,240,245,250,255,260,265,270,275,280,285,290,300]
cbar = plt.colorbar(img, cmap=cmap, norm=norm, boundaries=bounds, ticks=[200,205,210,215,220,225,230,235,240,245,250,255,260,265,270,275,280,285,290,300])
cbar.ax.set_yticklabels([200,205,210,215,220,225,230,235,240,245,250,255,260,265,270,275,280,285,290,300], fontsize=10)

plt.title('cloud top temperature 1km', fontsize=10)
plt.savefig("cloud_top_temperature_1km_no_projection.png", dpi=200)
#plt.show()

About / Au sujet de

Created:
19 janvier 2016 22:07:39

Updated:
19 janvier 2016 22:07:39

License / Licence

MIT License