#!/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-GEOPROF/'

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())

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

# Add file
f = addfile(‘D:/Temp/hdf/2010128055614\_21420\_CS\_2B-GEOPROF\_GRANULE\_P\_R04_E03.hdf‘)
# Read data
vname = ‘Radar_Reflectivity‘ v_data = f\[vname\]
data = v_data\[:,:\]
v_height = f\[‘Height‘\]
height = v_height\[0,:\]
time = f\[‘Profile_time‘\]\[:\]
lon = f\[‘Longitude‘\]\[:\]
lat = f\[‘Latitude‘\]\[:\] # Read attributes
long\_name = v\_data.attrvalue(‘long_name‘)\[0\]
scale_factor = v_data.attrvalue(‘factor‘)\[0\]
valid_min = v_data.attrvalue(‘valid_range‘)\[0\]
valid_max = v_data.attrvalue(‘valid_range‘)\[1\]
units = v_data.attrvalue(‘units‘)\[0\]
units_h = v_height.attrvalue(‘units‘)\[0\] # Apply scale factor
valid\_max = valid\_max / scale_factor
valid_min = valid_min / scale_factor
data = data / scale_factor
data\[data>valid_max\] = nan
data\[data<valid_min\] = nan
data = transpose(data)
data = data\[::-1,:\] # Make a split window plot
subplot(2, 1, 1) # Contour the data
levs = arange(-38, 50, 2)
layer = imshow(time, height\[::-1\], data, levs)
colorbar(layer)
title(‘Radar Reflectivity Factor‘)
xlabel(‘Seconds since the start of the granule. (seconds)‘)
ylabel(‘Height (m)‘) # The 2nd plot is the trajectory
subplot(2, 1, 2)
axesm()
lworld = shaperead(‘D:/Temp/map/country1.shp‘)
geoshow(lworld, edgecolor=‘k‘)
plotm(lon, lat, ‘-b‘, linewidth=4) #scatterm(lon, lat, lon, size=4, edge=False, facecolor=‘b‘)
scatterm(lon\[0\], lat\[0\], size=6, facecolor=‘r‘)
xlim(-180, 180)
ylim(-90, 90)
title(‘Trajectory of Flight Path (starting point in red)‘)