Notebook Source Codes

Share your source codes

Source codes: 131

Create modis caliop cloudsat matched files for a month with python 3

					
#!/usr/bin/env python

import os
import glob
import numpy as np
import calendar
import re

#----------------------------------------------------------------------------------------#
# User inputs

year = 2008
month = 7

machine = 'nasa' # icare

create_database_directory = False

#----------------------------------------------------------------------------------------#
# local machine (paths)

if machine == 'nasa':

	path_to_myd021km = '/Volumes/Data/MYD021KM/'
	path_to_myd03 = '/Volumes/Data/Myd03/'
	path_to_myd35 = '/Volumes/DATA_03/modis/Myd35/'
	path_to_myd06_c5 = '/Volumes/DATA_03/modis/Myd06_c5/'
	path_to_myd06_c6 = '/Volumes/DATA_03/modis/Myd06_c6/'
	path_to_01km_clay = '/Volumes/Data/Caliop/v4/01km_CLay/'
	path_to_05km_clay = '/Volumes/Data/Caliop/v4/05km_CLay/'
	path_to_cloudsat = '/Volumes/Data/Cloudsat/'
	path_to_2b_cldclass_lidar = '/Volumes/Data/2B-CLDCLASS-LIDAR/'

	dest_path = '/Users/marchant/Desktop/modis_based_match_files/'

#----------------------------------------------------------------------------------------#
# icare server (paths)

if machine == 'icare':

	path_to_myd021km = '/DATA/LIENS/MODIS/MYD021KM.006/'
	path_to_myd03 = '/DATA/LIENS/MODIS/MYD03/'
	path_to_myd35 = '/DATA/LIENS/MODIS/MYD35_L2.006/'
	path_to_myd06_c5 = '/DATA/LIENS/MODIS/MYD06_L2.005/'
	path_to_myd06_c6 = '/DATA/LIENS/MODIS/MYD06_L2.006/'
	path_to_01km_clay = '/DATA/LIENS/CALIOP/01kmCLay.v4.10/'
	path_to_05km_clay = '/DATA/LIENS/CALIOP/05kmCLay.v4.10/'
	path_to_cloudsat = '/DATA/LIENS/CLOUDSAT/2B-GEOPROF/'
	path_to_2b_cldclass_lidar = '/DATA/LIENS/CLOUDSAT/2B-CLDCLASS-LIDAR.v04.02/'

	dest_path = '/home/benjamin.hg.marchant/data/modis_based_match_files/'

#----------------------------------------------------------------------------------------#
# Create directory (Optional)

if create_database_directory:

	os.system( 'mkdir ' + dest_path + str(year) )

	for day in np.arange(calendar.monthrange(year,month)[1]):
		d = './' + str(year) + '/' + str(year) + '_' + "%02d" % month + \
			 '_' + "%02d" % (day+1) + '/'
		os.system('mkdir ' + dest_path + d)

#----------------------------------------------------------------------------------------#
# create file

day_of_year_0 = 0
for i in range(month-1):
	day_of_year_0 = day_of_year_0 + calendar.monthrange(year,i+1)[1]

for day in np.arange(calendar.monthrange(year,month)[1]):
#for day in [0]: #np.arange(calendar.monthrange(year,month)[1]):

	day_of_year = day_of_year_0 + day + 1

	#------------------------------------------------------------------------------------# 
	# create paths

	d = str(year) + '/' + str(year) + '_' + "%02d" % month + '_' + "%02d" % (day+1) + '/'

	path_to_myd021km_day = path_to_myd021km + d
	path_to_myd03_day = path_to_myd03 + d
	path_to_myd35_day = path_to_myd35 + d
	path_to_myd06_c5_day = path_to_myd06_c5 + d
	path_to_myd06_c6_day = path_to_myd06_c6 + d			
	path_to_01km_clay_day = path_to_01km_clay + d
	path_to_05km_clay_day = path_to_05km_clay + d
	path_to_cloudsat_day = path_to_cloudsat + d
	path_to_2b_cldclass_lidar_day = path_to_2b_cldclass_lidar + d

	#------------------------------------------------------------------------------------# 
	# get all file names

	list_of_myd021km_files = glob.glob(path_to_myd021km_day + '/*hdf' )
	list_of_myd03_files = glob.glob(path_to_myd03_day + '/*hdf' )
	list_of_myd35_files = glob.glob(path_to_myd35_day + '/*hdf' )
	list_of_myd06_c5_files = glob.glob(path_to_myd06_c5_day + '/*hdf' )
	list_of_myd06_c6_files = glob.glob(path_to_myd06_c6_day + '/*hdf' )
	list_of_01km_clay_files = glob.glob(path_to_01km_clay_day + '/*hdf' )
	list_of_05km_clay_files = glob.glob(path_to_05km_clay_day + '/*hdf' )	
	list_of_cloudsat_files = glob.glob(path_to_cloudsat_day + '/*hdf' ) 
	list_of_2b_cldclass_lidar_files = glob.glob(path_to_2b_cldclass_lidar_day + '/*hdf' )  

	#------------------------------------------------------------------------------------# 
	# create file

	myd021km_name = 'none'		
	myd03_name = 'none'	
	myd35_name = 'none'					
	myd06_c5_name = 'none'				
	myd06_c6_name = 'none'				
	cal_1km_file = 'none'	
	cal_5km_file = 'none'	
	cloudsat_2B_geoprof_file = 'none'	
	cloudsat_2B_cldclass_lidar_file = 'none'	

	for modis_granule_name in list_of_myd03_files:
		
		myd03_name = modis_granule_name.split('/')[-1]
		
		hm = myd03_name.split('.')[2]

		modis_granule_hour = int(hm[0:2])
		modis_granule_minute = int(hm[2:4])

		#--------------------------------------------------------------------------------#
		# Caliop 01km_clay

		min = 200

		file_exists = False

		for caliop_01km_granule_name in list_of_01km_clay_files:
			
			caliop_granule_name = caliop_01km_granule_name.split('/')[-1]
			
			caliop_re = re.compile(str(year) + '-' + "%02d" % month + '-' + \
			                       "%02d" % (day+1) + "T(.*)-(.*)-(.*).hdf")

			obj = caliop_re.search(caliop_granule_name)

			caliop_granule_hour = int(obj.group(1))
			caliop_granule_minute = int(obj.group(2))
			
			diff = int(modis_granule_hour) * 60 + int(modis_granule_minute) \
			        - ( int(caliop_granule_hour) * 60 + int(caliop_granule_minute) )

			if diff < min and diff > 0 :
				min = diff
				cal_1km_file = caliop_granule_name
				file_exists = True

		#--------------------------------------------------------------------------------#
		# Caliop 05km_clay
			
		min = 200

		file_exists = False

		for caliop_05km_granule_name in list_of_05km_clay_files:
			
			caliop_granule_name = caliop_05km_granule_name.split('/')[-1]
			
			caliop_re = re.compile(str(year) + '-' + "%02d" % month + '-' + \
			                       "%02d" % (day+1) + "T(.*)-(.*)-(.*).hdf")

			obj = caliop_re.search(caliop_granule_name)

			caliop_granule_hour = int(obj.group(1))
			caliop_granule_minute = int(obj.group(2))

			diff = int(modis_granule_hour) * 60 + int(modis_granule_minute) \
			        - ( int(caliop_granule_hour) * 60 + int(caliop_granule_minute) )

			if diff < min and diff > 0 :
				min = diff
				cal_5km_file = caliop_granule_name
				file_exists = True

		#--------------------------------------------------------------------------------#
		# MYD021KM
		
		file_name_prefix = 'MYD021KM.A'+ str(year) + \
		                   "%03d" % day_of_year + \
		                   '.' + "%02d" % modis_granule_hour + \
		                   "%02d" % modis_granule_minute
		found_file = 0
		for myd021km_granule_name in list_of_myd021km_files:
			offset = len(path_to_myd021km_day)
			if file_name_prefix == myd021km_granule_name[offset:offset+22]:
				found_file = 1
				myd021km_name = myd021km_granule_name[offset:]
		if found_file == 0: myd021km_name = 'NoFile'

		#--------------------------------------------------------------------------------#
		# MYD35

		file_name_prefix = 'MYD35_L2.A'+ str(year) + \
		                   "%03d" % day_of_year + \
		                   '.' + "%02d" % modis_granule_hour + \
		                   "%02d" % modis_granule_minute
		found_file = 0
		for myd35_granule_name in list_of_myd35_files:
			offset = len(path_to_myd35_day)
			if file_name_prefix == myd35_granule_name[offset:offset+22]:
				found_file = 1
				myd35_name = myd35_granule_name[offset:]
		if found_file == 0: myd35_name = 'NoFile'

		#--------------------------------------------------------------------------------#
		# MYD06_c5

		file_name_prefix = 'MYD06_L2.A'+ str(year) + \
		                   "%03d" % day_of_year + \
		                   '.' + "%02d" % modis_granule_hour + \
		                   "%02d" % modis_granule_minute
		found_file = 0		
		for myd06_c5_granule_name in list_of_myd06_c5_files:
			offset = len(path_to_myd06_c5_day)
			if file_name_prefix == myd06_c5_granule_name[offset:offset+22]:
				found_file = 1			
				myd06_c5_name = myd06_c5_granule_name[offset:]
		if found_file == 0: myd06_c5_name = 'NoFile'		

		#--------------------------------------------------------------------------------#
		# MYD06_c6

		file_name_prefix = 'MYD06_L2.A'+ str(year) + \
		                   "%03d" % day_of_year + \
		                   '.' + "%02d" % modis_granule_hour + \
		                   "%02d" % modis_granule_minute
		found_file = 0		
		for myd06_c6_granule_name in list_of_myd06_c6_files:
			offset = len(path_to_myd06_c6_day)
			if file_name_prefix == myd06_c6_granule_name[offset:offset+22]:
				found_file = 1			
				myd06_c6_name = myd06_c6_granule_name[offset:]
		if found_file == 0: myd06_c6_name = 'NoFile'	

		#--------------------------------------------------------------------------------#
		# cloudsat (GEOPROFF)

		find_cloudsat_match_file = False
		min_minute_diff = 99
		for cloudsat_file in list_of_cloudsat_files:
			
			cloudsat_file_name = cloudsat_file.split('/')[-1]
			cloudsat_granule_hour = int(cloudsat_file_name[7:9])
			cloudsat_granule_minute = int(cloudsat_file_name[9:11])
			cloudsat_granule_second = int(cloudsat_file_name[11:13] )
			
			minute_diff = modis_granule_hour * 60 + modis_granule_minute - \
					( cloudsat_granule_hour * 60 + cloudsat_granule_minute )
			
			if (minute_diff < min_minute_diff ) and ( minute_diff > 0 ):
					min_minute_diff = minute_diff
					cloudsat_2B_geoprof_file = cloudsat_file_name
					find_cloudsat_match_file = True

		#--------------------------------------------------------------------------------#
		# cloudsat (2B CLDCLASS-LIDAR)

		find_2b_cldclass_lidar_match_file = False
		min_minute_diff = 99
		for cloudsat_file in list_of_2b_cldclass_lidar_files:

			cloudsat_file_name = cloudsat_file.split('/')[-1]
			cloudsat_granule_hour = int(cloudsat_file_name[7:9])
			cloudsat_granule_minute = int(cloudsat_file_name[9:11])
			cloudsat_granule_second = int(cloudsat_file_name[11:13] )

			minute_diff = modis_granule_hour * 60 + modis_granule_minute - \
					( cloudsat_granule_hour * 60 + cloudsat_granule_minute )

			if (minute_diff < min_minute_diff ) and ( minute_diff > 0 ):
					min_minute_diff = minute_diff
					cloudsat_2B_cldclass_lidar_file = cloudsat_file_name
					find_2b_cldclass_lidar_match_file = True

		#--------------------------------------------------------------------------------#
		# save result

		filename = 'modis_linked_files_' + \
		           "%02d" % modis_granule_hour + "%02d" % modis_granule_minute
		
		f = open(filename, 'w')

		f.write(myd021km_name + ' \n')			
		f.write(myd03_name + ' \n')	
		f.write(myd35_name + ' \n')				
		f.write(myd06_c5_name + ' \n')			
		f.write(myd06_c6_name + ' \n')			
		f.write(cal_1km_file + ' \n')
		f.write(cal_5km_file + ' \n')
		f.write(cloudsat_2B_geoprof_file + ' \n')
		f.write(cloudsat_2B_cldclass_lidar_file + ' \n')
		
		f.closed
		
		command = 'mv ' + filename + ' ' + dest_path + d
		
		os.system(command)

Fortran 90 module to read CloudSat 2B CLDCLASS products

					
MODULE get_cldclass_lidar_data

implicit none

!----------------------------------------------------------------------------------------!
! Dynamic Allocation of Arrays

integer*2, dimension(:,:), allocatable :: cldclass_lidar_Height
real(kind=4), dimension(:,:), allocatable :: cldclass_lidar_CloudLayerBase
integer*1, dimension(:,:), allocatable :: cldclass_lidar_LayerBaseFlag
real(kind=4), dimension(:,:), allocatable :: cldclass_lidar_CloudLayerTop
integer*1, dimension(:,:), allocatable :: cldclass_lidar_LayerTopFlag
real(kind=4), dimension(:,:), allocatable :: cldclass_lidar_CloudFraction
integer*1, dimension(:,:), allocatable :: cldclass_lidar_CloudPhase
integer*1, dimension(:,:), allocatable :: cldclass_lidar_CloudPhaseConfidenceLevel
integer*1, dimension(:,:), allocatable :: cldclass_lidar_CloudLayerType
real(kind=4), dimension(:,:), allocatable :: cldclass_lidar_CloudTypeQuality
integer*1, dimension(:,:), allocatable :: cldclass_lidar_PrecipitationFlag
integer*1, dimension(:,:), allocatable :: cldclass_lidar_Phase_log
real(kind=4), dimension(:,:), allocatable :: cldclass_lidar_Water_layer_top

!----------------------------------------------------------------------------------------!
! User Inputs (user_selected_sds)

integer :: nb_user_selected_sds
parameter ( nb_user_selected_sds = 13 )

character*300 selected_sds( nb_user_selected_sds )      

data selected_sds  / 'Height', &
& 'CloudLayerBase', &
& 'LayerBaseFlag', &
& 'CloudLayerTop', &
& 'LayerTopFlag', &
& 'CloudFraction', &
& 'CloudPhase', &
& 'CloudPhaseConfidenceLevel', &
& 'CloudLayerType', &
& 'CloudTypeQuality', &
& 'PrecipitationFlag', &
& 'Phase_log', &
& 'Water_layer_top' /

!----------------------------------------------------------------------------------------!
! Subroutine read_cldclass_lidar

CONTAINS

subroutine read_cldclass_lidar(cldclass_lidar_file_name)

use read_write_hdf_file_variables

implicit none

integer :: i,j,k,l,m,n
integer :: AllocateStatus, DeAllocateStatus
integer :: ierr ! etat de reussite (si ==0) ou echec (si /= 0)
character*300 :: cldclass_lidar_file_name

integer sfstart, sfselect, sfrdata, sfendacc, sfend , sffattr , sfrattr, var_id, attr_id, err_code
integer status, n_attrs
integer sffinfo , sfginfo

character*300 vdata_name 
integer file_id, vdata_id, vdata_ref
integer n_records
character*80 fieldname_lists 
integer FULL_INTERLACE
parameter (FULL_INTERLACE = 0)
integer interlace_mode, vdata_size

!----------------------------------------------------------------------------------------!
! SDS

sd_id = sfstart(cldclass_lidar_file_name, DFACC_READ)
status = sffinfo( sd_id , n_datasets , n_file_attrs )

do index = 0, n_datasets - 1

	sds_id = sfselect(sd_id, index)
	status = sfginfo( sds_id , name , rank , dim_sizes , data_type , n_attrs )
    sds_name = name(1:len(name))

if ( ANY( selected_sds==sds_name ) ) then

	do i = 1 , rank
	start(i) = 0
	edges(i) = dim_sizes(i)
	stride(i) = 1
	end do

	if( sds_name == 'Height' )then
	allocate( cldclass_lidar_Height(dim_sizes(1),dim_sizes(2)), STAT = AllocateStatus)
	IF (AllocateStatus /= 0) STOP "*** Not enough memory (AllocateStatus cldclass_lidar_Height)***"
	status = sfrdata( sds_id, start, stride, edges, cldclass_lidar_Height )
	end if	
	
	if( sds_name == 'CloudLayerBase' )then
	allocate( cldclass_lidar_CloudLayerBase(dim_sizes(1),dim_sizes(2)), STAT = AllocateStatus)
	IF (AllocateStatus /= 0) STOP "*** Not enough memory (AllocateStatus cldclass_lidar_CloudLayerBase)***"
	status = sfrdata( sds_id, start, stride, edges, cldclass_lidar_CloudLayerBase )
	end if	
	
	if( sds_name == 'LayerBaseFlag' )then
	allocate( cldclass_lidar_LayerBaseFlag(dim_sizes(1),dim_sizes(2)), STAT = AllocateStatus)
	IF (AllocateStatus /= 0) STOP "*** Not enough memory (AllocateStatus cldclass_lidar_LayerBaseFlag)***"
	status = sfrdata( sds_id, start, stride, edges, cldclass_lidar_LayerBaseFlag )
	end if	
	
	if( sds_name == 'CloudLayerTop' )then
	allocate( cldclass_lidar_CloudLayerTop(dim_sizes(1),dim_sizes(2)), STAT = AllocateStatus)
	IF (AllocateStatus /= 0) STOP "*** Not enough memory (AllocateStatus cldclass_lidar_CloudLayerTop)***"
	status = sfrdata( sds_id, start, stride, edges, cldclass_lidar_CloudLayerTop )
	end if	
	
	if( sds_name == 'LayerTopFlag' )then
	allocate( cldclass_lidar_LayerTopFlag(dim_sizes(1),dim_sizes(2)), STAT = AllocateStatus)
	IF (AllocateStatus /= 0) STOP "*** Not enough memory (AllocateStatus cldclass_lidar_LayerTopFlag)***"
	status = sfrdata( sds_id, start, stride, edges, cldclass_lidar_LayerTopFlag )
	end if	
	
	if( sds_name == 'CloudFraction' )then
	allocate( cldclass_lidar_CloudFraction(dim_sizes(1),dim_sizes(2)), STAT = AllocateStatus)
	IF (AllocateStatus /= 0) STOP "*** Not enough memory (AllocateStatus cldclass_lidar_CloudFraction)***"
	status = sfrdata( sds_id, start, stride, edges, cldclass_lidar_CloudFraction )
	end if	
	
	if( sds_name == 'CloudPhase' )then
	allocate( cldclass_lidar_CloudPhase(dim_sizes(1),dim_sizes(2)), STAT = AllocateStatus)
	IF (AllocateStatus /= 0) STOP "*** Not enough memory (AllocateStatus cldclass_lidar_CloudPhase)***"
	status = sfrdata( sds_id, start, stride, edges, cldclass_lidar_CloudPhase )
	end if	
	
	if( sds_name == 'CloudPhaseConfidenceLevel' )then
	allocate( cldclass_lidar_CloudPhaseConfidenceLevel(dim_sizes(1),dim_sizes(2)), STAT = AllocateStatus)
	IF (AllocateStatus /= 0) STOP "*** Not enough memory (AllocateStatus cldclass_lidar_CloudPhaseConfidenceLevel)***"
	status = sfrdata( sds_id, start, stride, edges, cldclass_lidar_CloudPhaseConfidenceLevel )
	end if	
	
	if( sds_name == 'CloudLayerType' )then
	allocate( cldclass_lidar_CloudLayerType(dim_sizes(1),dim_sizes(2)), STAT = AllocateStatus)
	IF (AllocateStatus /= 0) STOP "*** Not enough memory (AllocateStatus cldclass_lidar_CloudLayerType)***"
	status = sfrdata( sds_id, start, stride, edges, cldclass_lidar_CloudLayerType )
	end if	
	
	if( sds_name == 'CloudTypeQuality' )then
	allocate( cldclass_lidar_CloudTypeQuality(dim_sizes(1),dim_sizes(2)), STAT = AllocateStatus)
	IF (AllocateStatus /= 0) STOP "*** Not enough memory (AllocateStatus cldclass_lidar_CloudTypeQuality)***"
	status = sfrdata( sds_id, start, stride, edges, cldclass_lidar_CloudTypeQuality )
	end if	
	
	if( sds_name == 'PrecipitationFlag' )then
	allocate( cldclass_lidar_PrecipitationFlag(dim_sizes(1),dim_sizes(2)), STAT = AllocateStatus)
	IF (AllocateStatus /= 0) STOP "*** Not enough memory (AllocateStatus cldclass_lidar_PrecipitationFlag)***"
	status = sfrdata( sds_id, start, stride, edges, cldclass_lidar_PrecipitationFlag )
	end if	
	
	if( sds_name == 'Phase_log' )then
	allocate( cldclass_lidar_Phase_log(dim_sizes(1),dim_sizes(2)), STAT = AllocateStatus)
	IF (AllocateStatus /= 0) STOP "*** Not enough memory (AllocateStatus cldclass_lidar_Phase_log)***"
	status = sfrdata( sds_id, start, stride, edges, cldclass_lidar_Phase_log )
	end if	
	
	if( sds_name == 'Water_layer_top' )then
	allocate( cldclass_lidar_Water_layer_top(dim_sizes(1),dim_sizes(2)), STAT = AllocateStatus)
	IF (AllocateStatus /= 0) STOP "*** Not enough memory (AllocateStatus cldclass_lidar_Water_layer_top)***"
	status = sfrdata( sds_id, start, stride, edges, cldclass_lidar_Water_layer_top )
	end if	

end if

end do

!----------------------------------------------------------------------------------------!

end subroutine read_cldclass_lidar

!----------------------------------------------------------------------------------------!

END MODULE get_cldclass_lidar_data

Read 2B-CLDCLASS-LIDAR SDS names and attributes using python 3 and pyhdf

					
#!/usr/bin/env python

from pyhdf.SD import SD, SDC 
from pyhdf.HDF import *
from pyhdf.VS import *

import pprint

#----------------------------------------------------------------------------------------#

file_path = './products/2B-CLDCLASS-LIDAR/'
file_name = '2008183012329_11573_CS_2B-CLDCLASS-LIDAR_GRANULE_P_R04_E02.hdf'

#----------------------------------------------------------------------------------------#
# Read HDF Files (VD data) 

f = HDF(file_path+file_name, SDC.READ) 
vs = f.vstart() 

print(vs.vdatainfo() )

Latitude = vs.attach('Latitude')
Longitude = vs.attach('Longitude')

a = Latitude[:]

print('VD data')
print(len(a),Latitude[0])

print( Latitude.attrinfo() )

print(min(a),max(a))

for i in range(20000,20020):
	print(a[i])

Latitude.detach() # "close" the vdata
Longitude.detach() # "close" the vdata
vs.end() # terminate the vdata interface
f.close() 

#----------------------------------------------------------------------------------------#

file = SD(file_path+file_name, SDC.READ)

file_info = file.info()
print(file_info)  # number of sds and metadata

#----------------------------------------------------------------------------------------#
# print sds names

datasets_dic = file.datasets()

#for idx,sds in enumerate(datasets_dic.keys()):
#	print idx,sds 

sds_dic = {}
for key, value in datasets_dic.items():
	#print( key, value, value[3] )
	sds_dic[value[3]] = key

pprint.pprint( sds_dic )

#----------------------------------------------------------------------------------------#

datasets_dic = file.datasets()

for key, value in datasets_dic.items():
	#print( key, value, value[3] )
	print('---------------------------')
	print('------------>', key)
	sds_obj = file.select(key) # select sds
	data = sds_obj.get() 
	sds_info = sds_obj.info()
	print( 'sds info' )
	print( sds_info )
	print( sds_info[0], sds_info[1] )
	print( 'sds attributes' )
	pprint.pprint( sds_obj.attributes() )

Create global maps from a Modis Caliop CloudSat sample using python

					
#!/usr/bin/env python

from mpl_toolkits.basemap import Basemap, cm

from scipy import stats, mgrid, c_, reshape, random, rot90
from pylab import *
from mpl_toolkits.axes_grid1 import make_axes_locatable

import numpy as np
import matplotlib.pyplot as plt
import seaborn as sn
import pandas as pd

import matplotlib as mpl

mpl.style.use('seaborn')

input_filename = '2008_01_sample_00.txt'

data = np.loadtxt(input_filename)

month_list = ['January','February ','March','April','May','June',
			  'July','August','September','October','November','December']

year = input_filename[0:4]
month = month_list[int(input_filename[5:7])-1]

#---------------------------------------------------------------------------------------#

data = data[~np.all(data == 0, axis=1)]

nb_tot_obs = len(data[:,0])

print('NUmber of observations ----> ', nb_tot_obs)

lat_min = np.min(data[:,0])
lat_max = np.max(data[:,0])

fig_count = 0 

#---------------------------------------------------------------------------------------#
# Scatter Plots Observations

fig = plt.figure()

ax = fig.add_subplot(111)

m = Basemap(projection='cyl',
			llcrnrlat=-90,urcrnrlat=90, llcrnrlon=-180,urcrnrlon=180,
			resolution='c')
m.drawcoastlines()
m.drawparallels(np.arange(-90.,90.,30.))
m.drawmeridians(np.arange(-180.,180.,30.))

m.scatter(data[:,1],data[:,0],s=0.02)

plt.title('Obs. Scatter Points \n Modis Caliop CloudSat Dataset (daytime only) \
		   \n '+ month + ' ' + year,fontsize=10)

plt.xlim(-180,180)
plt.ylim(-90,90)

fig_count = fig_count + 1
filename ='monthly_stats_sample_' + '{:02d}'.format(fig_count) + '.png'
plt.savefig(filename, dpi=200, bbox_inches='tight' )
#plt.show()
plt.close()

#---------------------------------------------------------------------------------------#
# Create kde Function

def create_kde(data_loc):
	xmin,xmax = -90.0, 90.0
	ymin,ymax = -180.0,180.0
	X, Y = np.mgrid[xmin:xmax:180j, ymin:ymax:360j]
	positions = np.vstack([X.ravel(), Y.ravel()])
	values = np.vstack([data_loc[:,0],data_loc[:,1]])
	kernel = stats.gaussian_kde(values)
	output_kde = np.reshape(kernel(positions).T, X.shape)
	return output_kde

#---------------------------------------------------------------------------------------#
# Create plot kde global map function

def create_kde_map(kde_img,fig_title,fig_count):
	fig = plt.figure()
	ax = fig.add_subplot(111)
	m = Basemap(projection='cyl',llcrnrlat=-90,urcrnrlat=90, \
				llcrnrlon=-180,urcrnrlon=180,resolution='c')
	m.drawcoastlines()
	m.drawparallels(np.arange(-90.,90.,30.))
	m.drawmeridians(np.arange(-180.,180.,30.))
	img= m.imshow(kde_img, interpolation='nearest', origin='lower', \
				  extent=[-180,180,-90,90], cmap=cm.jet)
	plt.title(fig_title,fontsize=10)
	divider = make_axes_locatable(ax)
	cax = divider.append_axes("right", size="5%", pad=0.1)
	cbar = fig.colorbar(img,format='%1.2g', cax=cax)
	cbar.set_label('')
	filename ='monthly_stats_sample_' + '{:02d}'.format(fig_count) + '.png'
	plt.savefig(filename, dpi=200, bbox_inches='tight' )
	plt.close()

#---------------------------------------------------------------------------------------#
# Create plot kde global map function

def create_fraction_map(fraction_img,fig_title,fig_count):
	fig = plt.figure()
	ax = fig.add_subplot(111)
	m = Basemap(projection='cyl',llcrnrlat=-90,urcrnrlat=90, \
				llcrnrlon=-180,urcrnrlon=180,resolution='c')
	m.drawcoastlines()
	m.drawparallels(np.arange(-90.,90.,30.))
	m.drawmeridians(np.arange(-180.,180.,30.))
	colors = [(0.75,0.75,0.75)] + [(cm.jet(i)) for i in range(1,256)]
	new_map = matplotlib.colors.LinearSegmentedColormap.from_list('new_map',colors, N=256)
	bounds = [-10.0,0.0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0]
	norm = matplotlib.colors.BoundaryNorm(bounds, new_map.N)
	img= m.imshow(fraction_img, interpolation='nearest', cmap=new_map, origin='lower', \
				  norm=norm,extent=[-180,180,-90,90],vmin=0,vmax=1)
	plt.title(fig_title,fontsize=10)
	divider = make_axes_locatable(ax)
	cax = divider.append_axes("right", size="5%", pad=0.1)
	ticks_at = [0.0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0]
	cbar = fig.colorbar(img,ticks=ticks_at,format='%1.2g', cax=cax)
	cbar.set_label('')
	filename ='monthly_stats_sample_' + '{:02d}'.format(fig_count) + '.png'
	plt.savefig(filename, dpi=200, bbox_inches='tight' )
	plt.close()

#---------------------------------------------------------------------------------------#
# KDE Observations

kde_obs = create_kde(data)

fig_count = fig_count + 1

fig_title = 'PDF (obs) \n Modis Caliop CloudSat Dataset (day only)\n'+ month + ' ' + year

create_kde_map(kde_obs,fig_title,fig_count)

#---------------------------------------------------------------------------------------#
# MODIS Cloud Fraction

data_modis_cloudy = data[(data[:,2] == 1) & (data[:,3] < 2)]

nb_tot_modis_cloudy = len(data_modis_cloudy[:,0])

kde_modis_cloudy = create_kde(data_modis_cloudy)

fraction_img = ( kde_modis_cloudy * nb_tot_modis_cloudy ) / ( kde_obs * nb_tot_obs )

for i in range(180):
	for j in range(360):
		if i < (90 + lat_min): fraction_img[i,j] = -9.0
		if i > (90 + lat_max): fraction_img[i,j] = -9.0

fig_count = fig_count + 1

fig_title = 'P( Modis Cloudy | $(\lambda,\phi)$)' + \
		   ' \n Modis Caliop CloudSat Dataset (day only) \n '+ month + ' ' + year

create_fraction_map(fraction_img,fig_title,fig_count)

#---------------------------------------------------------------------------------------#
# Caliop (1km & 5km) + CloudSat Cloud Fraction

data_cc_cloudy = data[data[:,5] > 0]

nb_tot_cc_cloudy= len(data_cc_cloudy[:,0])

kde_cc_cloudy = create_kde(data_cc_cloudy)

fraction_img = ( kde_cc_cloudy * nb_tot_cc_cloudy ) / ( kde_obs * nb_tot_obs )

for i in range(180):
	for j in range(360):
		if i < (90 + lat_min): fraction_img[i,j] = -9.0
		if i > (90 + lat_max): fraction_img[i,j] = -9.0

fig_count = fig_count + 1

fig_title = 'P( C-C Cloudy | $(\lambda,\phi)$)' + \
		   ' \n Modis Caliop CloudSat Dataset (day only) \n '+ month + ' ' + year

create_fraction_map(fraction_img,fig_title,fig_count)

#---------------------------------------------------------------------------------------#
# Modis + Caliop (1km & 5km) + CloudSat Cloud Fraction

data_filtered = data[(data[:,2] == 1) & (data[:,3] < 2) & (data[:,5] > 0)]

nb_tot_data_filtered = len(data_filtered[:,0])

kde_data_filtered = create_kde(data_filtered)

fraction_img = ( kde_data_filtered * nb_tot_data_filtered ) / ( kde_obs * nb_tot_obs )

for i in range(180):
	for j in range(360):
		if i < (90 + lat_min): fraction_img[i,j] = -9.0
		if i > (90 + lat_max): fraction_img[i,j] = -9.0

fig_count = fig_count + 1

fig_title = 'P( C-C & Modis Cloudy | $(\lambda,\phi)$)' + \
		   ' \n Modis Caliop CloudSat Dataset (day only) \n '+ month + ' ' + year

create_fraction_map(fraction_img,fig_title,fig_count)

#---------------------------------------------------------------------------------------#
# Caliop Cloudsat Clear vs MODIS cloudy Fraction

data_filtered = data[(data[:,2] == 1) & (data[:,3] < 2) & (data[:,5] == 0)]

nb_tot_data_filtered = len(data_filtered[:,0])

kde_data_filtered = create_kde(data_filtered)

fraction_img = ( kde_data_filtered * nb_tot_data_filtered ) / ( kde_obs * nb_tot_obs )

for i in range(180):
	for j in range(360):
		if i < (90 + lat_min): fraction_img[i,j] = -9.0
		if i > (90 + lat_max): fraction_img[i,j] = -9.0

fig_count = fig_count + 1

fig_title = 'P( C-C Clear & Modis Cloudy | $(\lambda,\phi)$)' + \
		   ' \n Modis Caliop CloudSat Dataset (day only) \n '+ month + ' ' + year

create_fraction_map(fraction_img,fig_title,fig_count)

#---------------------------------------------------------------------------------------#
# Caliop Cloudsat Cloudy vs MODIS Clear Fraction

data_filtered = data[(data[:,2] == 1) & (data[:,3] > 1) & (data[:,5] > 0)]

nb_tot_data_filtered = len(data_filtered[:,0])

kde_data_filtered = create_kde(data_filtered)

fraction_img = ( kde_data_filtered * nb_tot_data_filtered ) / ( kde_obs * nb_tot_obs )

for i in range(180):
	for j in range(360):
		if i < (90 + lat_min): fraction_img[i,j] = -9.0
		if i > (90 + lat_max): fraction_img[i,j] = -9.0

fig_count = fig_count + 1

fig_title = 'P( C-C Cloudy & Modis Clear | $(\lambda,\phi)$)' + \
		   ' \n Modis Caliop CloudSat Dataset (day only) \n '+ month + ' ' + year

create_fraction_map(fraction_img,fig_title,fig_count)

Download a month of CloudSat 2B-CLDCLASS-LIDAR products from icare

					
#!/usr/bin/env python

import numpy as np
import os
import calendar

from ftplib import FTP

#----------------------------------------------------------------------------------------#
# user inputs

year = 2008
month = 7

path_to_cloudsat = '/DATA/LIENS/CLOUDSAT/2B-CLDCLASS-LIDAR.v04.02/'

icare_username = 'john.doe'
icare_password = '********'

#----------------------------------------------------------------------------------------#
# Create directory 

os.system( 'mkdir ' + str(year) )

for day in np.arange(calendar.monthrange(year,month)[1]):
	d = './' + str(year) + '/' + str(year) + '_' + "%02d" % month + \
		 '_' + "%02d" % (day+1) + '/'
	os.system('mkdir ' + d)
	
#----------------------------------------------------------------------------------------#
# Opening ICARE FTP connection

print( 'Opening FTP connection' )

ftp = FTP('ftp.icare.univ-lille1.fr')
ftp.login(icare_username,icare_password)

for day in np.arange(calendar.monthrange(int(year),int(month))[1]):

	directory = path_to_cloudsat + str(year) + '/' + str(year) + '_' + \
				"%02d" % month + '_' + "%02d" % (day+1) + '/'	

	ftp.cwd(directory)

	yyy = []
	ftp.retrlines('NLST', yyy.append)
	#if len(yyy) != 288:
	#	print 'Probleme: donnees manquantes'
	for j in np.arange(len(yyy)):
		file_name = yyy[j]
		print(file_name)
		ftp.retrbinary('RETR ' + file_name, open(file_name, 'wb').write)
		command = 'mv ' + file_name + ' ' + './' + str(year) + '/' + str(year) + \
		'_' + "%02d" % month + '_' + "%02d" % (day+1) + '/' + file_name 
		os.system(command)

print( 'Closing FTP connection')
print( ftp.close())

Plot MODIS C6 Surface Pressure with python and matplotlib

					
#!/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 pprint
import calendar

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

file_name = './inputs/myd06_granules_viz/MYD06_L2.A2008015.1435.006.2013342100940.hdf'

#----------------------------------------------------------------------------------------#
# Extract Year, Month, Day, Hour and Minutes from file name

path = file_name.split('/')
name = path[len(path)-1]
list = name.split('.')

hhmm = list[2]
hh = int(hhmm[0:2])
mm = int(hhmm[2:4])

year = int(list[1][1:5])
dayofyear = int(list[1][5:8])

month = 1
day = dayofyear
while day - calendar.monthrange(year,month)[1] > 0 and month <= 12:
	day = day - calendar.monthrange(year,month)[1]
	month = month + 1

#print year, month, day, hh, mm

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

file = SD(file_name, SDC.READ)

data_selected_id = file.select('Surface_Pressure')

data = data_selected_id.get()

data_selected_attributes = data_selected_id.attributes()

scale_factor = data_selected_attributes['scale_factor']
add_offset = data_selected_attributes['add_offset']
fillvalue = data_selected_attributes['_FillValue']

#data[data == fillvalue] = 0

data_shape = data.shape

#----------------------------------------------------------------------------------------#
# Plot data

data = (data - add_offset) * scale_factor
data[data == (fillvalue - add_offset) * scale_factor] = 0

print(data.min(),data.max())

cmap = [(0.0,0.0,0.0)] + [(cm.jet(i)) for i in range(1,256)] 
cmap = mpl.colors.ListedColormap(cmap)

bounds = [0] + [i for i in np.linspace(data.min(),data.max(),10)]

norm = mpl.colors.BoundaryNorm(bounds, cmap.N)

img = plt.imshow(np.fliplr(data), cmap=cmap, norm=norm,
                  interpolation='none', origin='lower')

cbar_bounds = [0] + [i for i in np.linspace(data.min(),data.max(),10)]
cbar_ticks = [0] + [i for i in np.linspace(data.min(),data.max(),10)]               
               
cbar = plt.colorbar(img, cmap=cmap, norm=norm, boundaries=cbar_bounds, ticks=cbar_ticks)
cbar.ax.set_yticklabels([round(i,2) for i in cbar_ticks], fontsize=10)

l = [int(i) for i in np.linspace(0,data_shape[1],6)]
plt.xticks(l, [i for i in reversed(l)], rotation=0, fontsize=8 )

l = [int(i) for i in np.linspace(0,data_shape[0],9)]
plt.yticks(l, l, rotation=0, fontsize=8 )

title = 'MODIS Surface Pressure \n MYD06 C6 ({}-{:02d}-{:02d}; {:02d}h{:02d})'
plt.title(title.format( year, month, day, hh, mm), fontsize=10)

plt.savefig("myd06_c6_surface_pressure_data.png", 
            bbox_inches='tight', dpi=200)
plt.show()

Plot MODIS C6 Surface Temperature with python and matplotlib

					
#!/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 pprint
import calendar

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

file_name = './inputs/myd06_granules_viz/MYD06_L2.A2008015.1435.006.2013342100940.hdf'

#----------------------------------------------------------------------------------------#
# Extract Year, Month, Day, Hour and Minutes from file name

path = file_name.split('/')
name = path[len(path)-1]
list = name.split('.')

hhmm = list[2]
hh = int(hhmm[0:2])
mm = int(hhmm[2:4])

year = int(list[1][1:5])
dayofyear = int(list[1][5:8])

month = 1
day = dayofyear
while day - calendar.monthrange(year,month)[1] > 0 and month <= 12:
	day = day - calendar.monthrange(year,month)[1]
	month = month + 1

#print year, month, day, hh, mm

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

file = SD(file_name, SDC.READ)

data_selected_id = file.select('Surface_Temperature')

data = data_selected_id.get()

data_selected_attributes = data_selected_id.attributes()

scale_factor = data_selected_attributes['scale_factor']
add_offset = data_selected_attributes['add_offset']
fillvalue = data_selected_attributes['_FillValue']

#data[data == fillvalue] = 0

data_shape = data.shape

#----------------------------------------------------------------------------------------#
# Plot data

data = (data - add_offset) * scale_factor
data[data == (fillvalue - add_offset) * scale_factor] = 0

print(data.min(),data.max())

cmap = [(0.0,0.0,0.0)] + [(cm.jet(i)) for i in range(1,256)] 
cmap = mpl.colors.ListedColormap(cmap)

bounds = [0] + [i for i in np.linspace(data.min(),data.max(),10)]

norm = mpl.colors.BoundaryNorm(bounds, cmap.N)

img = plt.imshow(np.fliplr(data), cmap=cmap, norm=norm,
                  interpolation='none', origin='lower')

cbar_bounds = [0] + [i for i in np.linspace(data.min(),data.max(),10)]
cbar_ticks = [0] + [i for i in np.linspace(data.min(),data.max(),10)]               
               
cbar = plt.colorbar(img, cmap=cmap, norm=norm, boundaries=cbar_bounds, ticks=cbar_ticks)
cbar.ax.set_yticklabels([round(i,2) for i in cbar_ticks], fontsize=10)

l = [int(i) for i in np.linspace(0,data_shape[1],6)]
plt.xticks(l, [i for i in reversed(l)], rotation=0, fontsize=8 )

l = [int(i) for i in np.linspace(0,data_shape[0],9)]
plt.yticks(l, l, rotation=0, fontsize=8 )

title = 'MODIS Surface Temperature \n MYD06 C6 ({}-{:02d}-{:02d}; {:02d}h{:02d})'
plt.title(title.format( year, month, day, hh, mm), fontsize=10)

plt.savefig("myd06_c6_surface_temperature_data.png", 
            bbox_inches='tight', dpi=200)
plt.show()

How to plot MODIS L2 cloud phase infrared 5km with python 3

					
#!/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 calendar

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

file_name = './inputs/myd06_granules_viz/MYD06_L2.A2008015.1435.006.2013342100940.hdf'

#----------------------------------------------------------------------------------------#
# Extract Year, Month, Day, Hour and Minutes from file name

path = file_name.split('/')
name = path[len(path)-1]
list = name.split('.')

hhmm = list[2]
hh = int(hhmm[0:2])
mm = int(hhmm[2:4])

year = int(list[1][1:5])
dayofyear = int(list[1][5:8])

month = 1
day = dayofyear
while day - calendar.monthrange(year,month)[1] > 0 and month <= 12:
	day = day - calendar.monthrange(year,month)[1]
	month = month + 1

#print year, month, day, hh, mm

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

file = SD(file_name, SDC.READ)

data_selected_id = file.select('Cloud_Phase_Infrared')

data = data_selected_id.get()

data_shape = data.shape

#----------------------------------------------------------------------------------------#
# Plot data

cmap = [(0.0,0.0,0.0)] + [(cm.jet(i)) for i in range(1,256)] 
cmap = mpl.colors.ListedColormap(cmap)

bounds = [-0.5,0.5,1.5,2.5,3.5,4.5,5.5,6.5]

norm = mpl.colors.BoundaryNorm(bounds, cmap.N)

img = plt.imshow(np.fliplr(data), cmap=cmap, norm=norm,
                  interpolation='none', origin='lower')

cbar_bounds =  [-0.5,0.5,1.5,2.5,3.5,4.5,5.5,6.5]
cbar_ticks =  [0,1,2,3,4,5,6]               
cbar = plt.colorbar(img, cmap=cmap, norm=norm, boundaries=cbar_bounds, ticks=cbar_ticks)

yticklabels_list = ['0 - cloud free','1 - water cloud ','2 - ice cloud','3 - mixed cloud','4','5','6 - undetermined']
cbar.ax.set_yticklabels(yticklabels_list, fontsize=8)

l = [int(i) for i in np.linspace(0,data_shape[1],6)]
plt.xticks(l, [i for i in reversed(l)], rotation=0, fontsize=8 )

l = [int(i) for i in np.linspace(0,data_shape[0],9)]
plt.yticks(l, l, rotation=0, fontsize=8 )

title = 'Cloud Phase Infrared 5km \n MYD06 C6 ({}-{:02d}-{:02d}; {:02d}h{:02d})'

plt.title(title.format( year, month, day, hh, mm), fontsize=10)

plt.xticks(fontsize=8)
plt.yticks(fontsize=8)

plt.savefig("myd06_c6_cloud_phase_infrared.png", 
            bbox_inches='tight', dpi=200)
plt.show()

How to plot MODIS L2 cloud top height with python 3

					
#!/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 pprint
import calendar

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

file_name = './inputs/myd06_granules_viz/MYD06_L2.A2005099.1050.006.2014027131121.hdf'

#----------------------------------------------------------------------------------------#
# Extract Year, Month, Day, Hour and Minutes from file name

path = file_name.split('/')
name = path[len(path)-1]
list = name.split('.')

hhmm = list[2]
hh = int(hhmm[0:2])
mm = int(hhmm[2:4])

year = int(list[1][1:5])
dayofyear = int(list[1][5:8])

month = 1
day = dayofyear
while day - calendar.monthrange(year,month)[1] > 0 and month <= 12:
	day = day - calendar.monthrange(year,month)[1]
	month = month + 1

#print year, month, day, hh, mm

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

file = SD(file_name, SDC.READ)

data_selected_id = file.select('Cloud_Top_Height')

data = data_selected_id.get()

data_selected_attributes = data_selected_id.attributes()

fillvalue = data_selected_attributes['_FillValue']

data[data == fillvalue] = 0

data = data * 0.001

#----------------------------------------------------------------------------------------#
# Plot data

cmap = [(0.0,0.0,0.0)] + [(cm.jet(i)) for i in range(1,256)] 
cmap = mpl.colors.ListedColormap(cmap)

bounds = [0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18]

norm = mpl.colors.BoundaryNorm(bounds, cmap.N)

img = plt.imshow(np.fliplr(data), cmap=cmap, norm=norm,
                  interpolation='none', origin='lower')

cbar_bounds = [0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18]
cbar_ticks = [0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18]             
               
cbar = plt.colorbar(img, cmap=cmap, norm=norm, boundaries=cbar_bounds, ticks=cbar_ticks)
cbar.ax.set_yticklabels(cbar_ticks, fontsize=10)

plt.title('Cloud Top Height (Km) \n MYD06 C6 (2008-01-15; 14h35)', fontsize=10)

plt.savefig("myd06_c6_cloud_top_height.png", 
            bbox_inches='tight', dpi=200)
plt.show()

How to plot MODIS L2 cloud height method with python 3

					
#!/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 calendar

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

file_name = './inputs/myd06_granules_viz/MYD06_L2.A2005099.1050.006.2014027131121.hdf'

#----------------------------------------------------------------------------------------#
# Extract Year, Month, Day, Hour and Minutes from file name

path = file_name.split('/')
name = path[len(path)-1]
list = name.split('.')

hhmm = list[2]
hh = int(hhmm[0:2])
mm = int(hhmm[2:4])

year = int(list[1][1:5])
dayofyear = int(list[1][5:8])

month = 1
day = dayofyear
while day - calendar.monthrange(year,month)[1] > 0 and month <= 12:
	day = day - calendar.monthrange(year,month)[1]
	month = month + 1

#print year, month, day, hh, mm

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

file = SD(file_name, SDC.READ)

data_selected_id = file.select('Cloud_Height_Method')

data = data_selected_id.get()

data_selected_attributes = data_selected_id.attributes()

fillvalue = data_selected_attributes['_FillValue']

data[data == fillvalue] = 0

#----------------------------------------------------------------------------------------#
# Plot data

cmap = [(0.0,0.0,0.0)] + [(cm.jet(i)) for i in range(1,256)] 
cmap = mpl.colors.ListedColormap(cmap)

bounds = [-0.5,0.5,1.5,2.5,3.5,4.5,5.5,6.5]

norm = mpl.colors.BoundaryNorm(bounds, cmap.N)

img = plt.imshow(np.fliplr(data), cmap=cmap, norm=norm,
                  interpolation='none', origin='lower')

cbar_bounds =  [-0.5,0.5,1.5,2.5,3.5,4.5,5.5,6.5]
cbar_ticks =  [0,1,2,3,4,5,6]               
               
cbar = plt.colorbar(img, cmap=cmap, norm=norm, boundaries=cbar_bounds, ticks=cbar_ticks)
cbar.ax.set_yticklabels(cbar_ticks, fontsize=8)

title = 'Cloud Height Method \n MYD06 C6 ({}-{:02d}-{:02d}; {:02d}h{:02d})'
plt.title(title.format( year, month, day, hh, mm), fontsize=10)

plt.xticks(fontsize=8)
plt.yticks(fontsize=8)

plt.savefig("myd06_c6_cloud_height_method.png", 
            bbox_inches='tight', dpi=200)
plt.show()

About


To add a source code to your notebook just click on the link below. Sharing scientific source code allows to understand and reproduce your results:

Add a code

Search


Source codes can be retrieved using the following form: