#!/usr/bin/env python

import numpy as np
import numpy.matlib
import matplotlib.pyplot as plt
import matplotlib.cm as cm
import matplotlib.colors as mcolors
from netCDF4 import Dataset
import os as os
import resource
from mpl_toolkits.basemap import Basemap
#from scipy.stats import norm



def listevar(fichier):
	datago = Dataset(fichier, 'r')
	x = datago.variables.keys()
	print x

def extraction(fichier):
	datago = Dataset(fichier, 'r')
	lg = datago.variables['longitude'][:]
	lt = datago.variables['latitude'][:]
	ti = datago.variables['time'][:]
	atb = datago.variables['ATB'][:,:]
	SR = datago.variables['instant_SR'][:,:]
	sol_alt = datago.variables['SE'][:]
	alt = datago.variables['alt_bound'][:,:]
	datago.close()
	u, v = np.shape(alt)
	alt_bas = alt[0,:]
	alt_haut = alt[1,:]
	alt = np.zeros((v,1))
	for i in range(0,v):
		alt[i] = alt_bas[i]
	#alt[v] = alt_haut[v-1]
	alt = alt[:,0]
	return lg, lt, ti, atb,alt, sol_alt, SR



def affiche(n1,n2):
	lt = lt[n1:n2]
	sol_alt = sol_alt[n1:n2]
	u,v = np.shape(atb[n1:n2,:])
	trait = lt[int (((n2-n1)/2))] 
	plt.plot ([trait,trait] , [alt[0],alt[-1]],'w--',linewidth=2)
	plt.plot(lt, sol_alt,'w',linewidth=2)
	plt.pcolormesh(lt,alt,atb[n1:n2,:].T)
	plt.ylim(-1, 16)
	plt.xlim(lt[0],lt[-1])
	#plt.gca().invert_xaxis()
	plt.clim(0, 0.005)
	plt.xlabel('Latitude')
	plt.ylabel('Altitude [km]')
	plt.colorbar()	  
	plt.show() 

def afficheSR(lg,lt,ti,atb,alt,sol_alt):
	n1 = 0
	n2,x = np.shape(atb)
	lt = lt[n1:n2]
	sol_alt = sol_alt[n1:n2]
	u,v = np.shape(atb[n1:n2,:])	
 	sol_alt = sol_alt[n1:n2]
	fig = plt.figure(figsize=(9, 6), dpi=100)
	plt.plot(lt, sol_alt,'w',linewidth=2)
	plt.pcolormesh(lt,alt,atb[n1:n2,:].T)
	plt.gca().invert_xaxis()
	plt.ylim(-1, 16)
	plt.clim(1e-5, 5)
	plt.xlabel('Latitude')
	plt.ylabel('Altitude [km]')
	plt.colorbar()	  
	plt.show() 

def lonlat(n1, n2,lg,lt,ti,atb):
	from mpl_toolkits.basemap import Basemap
	m = Basemap()
	m.plot(lg, lt)
	lg = lg[n1:n2]
	lt = lt[n1:n2]
	m.plot(lg, lt, 'r')
	m.drawcoastlines()
	m.fillcontinents()
	plt.show()

def draw():
	m = Basemap()
	lon1 ,lon2,lat1,lat2 = -110,-130,-20,-40 # shallowcu
	#lon1 ,lon2,lat1,lat2 = 180,70,-15,15 # warmpool
	fig = plt.figure(figsize=(14, 7), dpi=100)
	m.plot([lon1,lon1], [lat1,lat2], lw = 2, color = 'Gold')
	m.plot([lon2,lon2], [lat1,lat2], lw = 2, color = 'Gold')
	m.plot([lon1,lon2], [lat1,lat1], lw = 2, color = 'Gold')
	m.plot([lon1,lon2], [lat2,lat2], lw = 2, color = 'Gold')
	m.bluemarble()
	parallels = np.arange(-90.,90,20.)
	m.drawparallels(parallels,labels=[1,0,0,0],fontsize=10,color = 'Chocolate')
	meridians = np.arange(-180.,180.,20.)
	m.drawmeridians(meridians,labels=[0,0,0,1],fontsize=10, color = 'Chocolate')
	#m.drawcoastlines()
	#m.fillcontinents()
	#titre = 'Warmpool'
	titre = 'Shallowcu'
	path = titre + '_carte'
	plt.title(titre)
	plt.savefig(path)
	#plt.show()
	plt.close()

def profil(n,lg,lt,ti,atb,alt,fichier):
	datago = Dataset(fichier, 'r')
	atmol = datago.variables['ATB_mol'][:]
	p1, = plt.semilogx(atb[n,:], alt)
	p2, = plt.semilogx(atmol[n,:], alt)
	plt.legend([p1,p2], ["ATB mesure","ATB moleculaire"])
	#plt.xlim(1e-7, 1e-1)
	plt.ylim(0,20)
	plt.xlabel('ATB [km-1 sr-1]')
	plt.ylabel('Altitude [km]')
	plt.title('Profil ' + str(lt[n]))
	plt.show()

def profil_SR(n,lg,lt,ti,atb,alt,fichier):
	plt.plot(atb[n,:], alt)
	p1, = plt.plot([5,5],[alt[0],alt[-1]],'--')
	p2, = plt.plot([0.01,0.01],[alt[0],alt[-1]], '-*')
	plt.legend([p1,p2], ["SR = 5","SR = 0,01"])
	#axes.set_xticks(ax)
	#plt.xscale('symlog', linthreshx=0.1)
	plt.xlim(-0.1, 50)
	plt.ylim(0,20)
	plt.xlabel('SR')
	plt.ylabel('Altitude [km]')
	plt.title('Profil ' + str(lt[n]))
	plt.show()
	#path = fichier
	#path += '.png'
	#plt.savefig(path)
	#plt.close()
	

def coord_pole(lg,lt,ti,atb, sol_alt):
	idx = ( lt < -60 ) & ( lt >= -90 )
	lt6090 = lt[idx]
	lg6090 = lg[idx]
	ti6090 = ti[idx]
	sol_alt = sol_alt[idx]
	atb = atb [idx,:]
	n = len(lt6090)
	return n, lt6090, lg6090, ti6090,atb, sol_alt

def alt_max(SR, alt, n):
	u, v = np.shape(SR)
	j = - 1
	if (SR[n,0] > 0.01 or SR[n,0] < -10) or (SR[n,1] > 0.01 or SR[n,1] < -10):
		return j
	i = 2
	while i < v:
		if (SR[n,i] > -10) and (SR[n,i] < 0.01):
			j = i
		elif (SR[n,i] > 0.01) and (SR[n,i-1] > 0.01):
			return j
		i = i + 1
	return j

def alt_max_orbite(alt, SR):
	p,v = np.shape(SR)
	tab1 = np.zeros((p))
	tab2 = np.zeros((p))
	for o in range (0,p):
		j = alt_max(SR,alt,o)
		if j != -1 :
			tab1[o] = alt[j]
			#if tab1[o] > 5.: print o
		else : tab1[o] = -10
		tab2[o] = o
	return tab1
	'''
	tab1 = np.ma.masked_where(tab1 <= -1, tab1)
	fig = plt.figure(figsize=(15, 10), dpi=300)
	plt.plot(lt,tab1)
	plt.plot(lt, sol_alt,'g',linewidth=2)
	plt.xlim(np.min(lt), np.max(lt))
	plt.ylabel('Altitude')
	plt.xlabel('indice')
	plt.show()
	'''
def alt_max_orbite_flag(alt, SR):
	u,v= np.shape(SR)
	tab1 = np.zeros((u,v))
	for o in range (0,u):
		j = alt_max(SR,alt,o)
		if j != -1 :tab1[:][o] = alt[j]
		else : tab1[:][o] = -10
	return tab1

def alt_min_SR5(alt_array,SR,altmax):
	u,v = np.shape(SR)
	alt_SR = np.zeros((u,v))
	idx = (SR > 5)
	alt_copy = alt_array.copy()
	
	alt_copy[~idx] = 9999.
	alt_copy[altmax < 0] = 9999.
	alt_SR = np.min(alt_copy, axis=1)
	alt_SR[alt_SR == 9999] = -1.
	return alt_SR

def alt_min_SR5_orbite(alt_array,SR,altmax):
	u,v = np.shape(SR)
	alt_SR = np.zeros((u,v))
	alt_tab = np.zeros((u,v))
	alt_total = np.zeros((u,v))
	idx = (SR > 5)
	alt_copy = alt_array.copy()
	
	alt_copy[~idx] = 9999.
	alt_copy[altmax < 0] = 9999.
	alt_SR = np.min(alt_copy, axis=1)
	alt_SR[alt_SR == 9999] = -1.
	
	for i in range(0,v-1):
		if i == 0: alt_tab = alt_SR
		alt_tab = np.vstack((alt_tab,alt_SR))	

	alt_tab = alt_tab.T
	print alt_tab
	print np.shape(alt_tab)
	return alt_tab
	
def alt_max_profil(alt,SR,n,fichier):
	j = alt_max(SR,alt,n)
	s = np.shape(alt)
	if j != -1 :alt_ext= alt[j]
	else : alt_ext= -10
	print alt_ext
	profil_SR(n,lg,lt,ti,SR,alt,fichier)
	
def affiche_flag(lg,lt,atb,alt,sol_alt):
	n1 = 0
	n2, x = np.shape(atb)
	lt = lt[n1:n2]
	sol_alt = sol_alt[n1:n2]
	u,v = np.shape(atb[n1:n2,:])
	fig = plt.figure(figsize=(14, 7), dpi=100)
	plt.plot(lt, sol_alt,'Maroon', alpha=1)
	cmap = cm.colors.ListedColormap(['White', 'Navy', 'Lime', 'Chocolate','Gold'])
	plot = plt.pcolormesh(lt,alt,atb.T,cmap=cmap)
	plt.ylim(0, 183)
	plt.xlim(-80, 80)
	plt.clim(0, 4)
	plt.xlabel('Latitude')
	plt.ylabel('Altitude en km')

	cb =plt.colorbar(ticks=[0, 1.3, 2, 2.75,3.5])
	cb.ax.set_yticklabels(['Clear Sky','Nuages','Faux eteint','Zone de transition','Vrai eteint'])  	  
	path = 'SR_orbite_vrai_faux'
	plt.savefig(path)
	plt.close()

def affiche_flag_grl(lg,lt,atb,alt,sol_alt):
	n1 = 0
	n2, x = np.shape(atb)
	lt = lt[n1:n2]
	sol_alt = sol_alt[n1:n2]
	u,v = np.shape(atb[n1:n2,:])
	fig = plt.figure(figsize=(14, 7), dpi=100)
	plt.plot(lt, sol_alt,'Maroon', alpha=1)
	cmap = cm.colors.ListedColormap(['White', 'Navy','Gold', 'Chocolate'])
	plot = plt.pcolormesh(lt,alt,atb.T,cmap=cmap)
	plt.ylim(0, 18)
	plt.xlim(-80, 80)
	plt.clim(0, 3)
	plt.xlabel('Latitude')
	plt.ylabel('Altitude en km')

	cb =plt.colorbar(ticks=[0, 1.3, 2,3])
	cb.ax.set_yticklabels(['Clear Sky','Nuages','Eteint','Transition'])  	  
	path = 'SR_orbite_GRL'
	plt.savefig(path)
	plt.close()

def flag_orbite_grl(alt,SR):
	u,v = np.shape(SR)
	flag = np.zeros((u,v))
	alt_array = np.zeros((u,v))
	alt_SR = np.zeros((u,v))
	alt_array[:][:] = alt
	alt_max = alt_max_orbite_flag(alt,SR)
	alt_SR = np.matlib.repmat(alt_min_SR5(alt_array,SR,alt_max),v,1).T
	print np.shape(alt_SR)
	#Nuage
	idx = (SR > 5)
	flag[idx] = 1

	#Vrai eteint
	idx4 = (SR < 0.01) & (SR> -10) 
	flag[idx4] = 2
	
	#Zone de transition
	idx3 =  (alt_array > alt_max) & (alt_array < alt_SR)
	flag[idx3] = 3

	return flag

def flag_orbite(alt,SR):
	u,v = np.shape(SR)
	flag = np.zeros((u,v))
	alt_array = np.zeros((u,v))
	alt_SR = np.zeros((u,v))
	alt_array[:][:] = alt
	alt_max = alt_max_orbite_flag(alt,SR)
	alt_SR = np.matlib.repmat(alt_min_SR5(alt_array,SR,alt_max),v,1).T
	print np.shape(alt_SR)
	#Nuage
	idx = (SR > 5)
	flag[idx] = 1

	#Vrai eteint
	idx4 = (alt_array < alt_max)&(SR < 0.01) & (SR> -10)
	flag[idx4] = 4

	#Faux eteint
	idx2 = ((alt_array > alt_max)  & (alt_max!= -10)&(SR < 0.01) & (SR> -10)) | ((alt_max == -10) & (SR < 0.01) & (SR> -10))
	flag[idx2] = 2
	
	#Zone de transition
	idx3 =  (alt_array > alt_max) & (alt_array < alt_SR)
	flag[idx3] = 3

	return flag
	
def lecture_dossier(repertoire):
	u = np.max(np.shape(os.listdir(repertoire)))
	
	return u, os.listdir(repertoire)

def occur(flag,alt):
	clouds = (flag == 1).sum(axis=0)
	faux_eteint = (flag == 2).sum(axis=0) 
	cloud_opa = (flag == 3).sum(axis=0)
	vrai_eteint = (flag == 4).sum(axis=0)
	return clouds, faux_eteint, cloud_opa, vrai_eteint

def freq_occur(flag,alt):
	u,v = np.shape(flag)
	clouds, faux_eteint, cloud_opa, vrai_eteint = occur(flag,alt)
	clouds, faux_eteint, cloud_opa, vrai_eteint = (1.)*clouds / u, (1.)*faux_eteint / u, (1.)*cloud_opa / u, (1.)*vrai_eteint / u
	fig = plt.figure(figsize=(9, 15), dpi=100)
	p1, = plt.plot(vrai_eteint,alt)
	sum1 = vrai_eteint + cloud_opa
	p2, = plt.plot(sum1,alt)
	sum2 = vrai_eteint + cloud_opa + clouds
	p3, = plt.plot(sum2,alt)
	sum3 = vrai_eteint + cloud_opa + clouds + faux_eteint
	p4, = plt.plot(sum3,alt)
	#plt.xlim(0.0001, 0.1)
	plt.legend([p1,p2,p3,p4], ["Nuages","Faux eteint","Zone de transition","Vrai eteint"], loc=1)
	plt.xlabel("Frequence d'occurence")
	plt.ylabel("Altitude en  km")
	plt.show()
	
def occur_in_region(flag,alt,lt,latrange):
	idx = (lt < latrange[1]) & (lt > latrange[0])
	flag = flag[idx]
	clouds = (flag == 1).sum(axis=0)
	faux_eteint = (flag == 2).sum(axis=0) 
	cloud_opa = (flag == 3).sum(axis=0)
	vrai_eteint = (flag == 4).sum(axis=0)
	return clouds, faux_eteint, cloud_opa, vrai_eteint

def freq_occur_region(dossier,test):
	nbr_fichier,fichier = lecture_dossier(dossier)	

	dic_opa_cloud = dict()
	dic_cloud = dict()
	dic_vrai = dict()
	dic_faux = dict()
	#nbr_fichier=1
	nbr_fichier=np.max(nbr_fichier)
	latbands = {'tropique' :[-20,20],'pole_nord' :[60,90],'pole_sud' : [-90,-60],'mid_sud' :[-50,-30],'mid_nord' :[30, 50 ]}

	if test == 1 :

		for i in range(0,nbr_fichier):	
			fichier3 = dossier + fichier[i]
			lg, lt, ti, atb,alt, sol_alt, SR = extraction(fichier3)
			flag = flag_orbite(alt,SR)
			u,v = np.shape(flag)
			print fichier[i]
			for region, tab_reg  in latbands.items():
				#flag = flag[tab_reg,:]
				clouds, faux_eteint, cloud_opa, vrai_eteint = occur_in_region(flag, alt, lt, latbands[region])
				if region in dic_opa_cloud:
					dic_opa_cloud[region] += 1.* cloud_opa/ (u * nbr_fichier)
					dic_cloud[region] += 1.* clouds/ (u * nbr_fichier)
					dic_vrai[region] += 1.* vrai_eteint/ (u * nbr_fichier)
					dic_faux[region] += 1.* faux_eteint/ (u * nbr_fichier)
				else:
					dic_opa_cloud[region] = 1.* cloud_opa/ (u * nbr_fichier)
					dic_cloud[region] = 1.* clouds/ (u * nbr_fichier)
					dic_vrai[region] = 1.* vrai_eteint/ (u * nbr_fichier)
					dic_faux[region] = 1.* faux_eteint/ (u * nbr_fichier)
		for	region in latbands.keys():
			output = region
			np.savez(output, dic_opa_cloud[region], dic_cloud[region],dic_vrai[region],dic_faux[region],alt)

	if test == 0 :

		for	region in latbands.keys():
			fil = region +'.npz'
			npz = np.load(fil)
			dic_opa_cloud = npz['arr_0']
			dic_cloud = npz['arr_1']
			dic_vrai = npz['arr_2']
			dic_faux = npz['arr_3']
			alt = npz['arr_4']
			affiche_freq_separee(alt,dic_opa_cloud,dic_cloud,dic_vrai,dic_faux,region,fichier,nbr_fichier)

	def affiche_freq(alt,op,cl,vr,fa,region,fichier,nbr_fichier):
		fig = plt.figure(figsize=(9, 15), dpi=100)
		p1, = plt.plot(op,alt)
		sum1 = vr + op
		p2, = plt.plot(sum1,alt)
		sum2 = sum1 + fa
		p3, = plt.plot(sum2,alt)
		sum3 = sum2 + cl
		p4, = plt.plot(sum3,alt)
		plt.legend([p1], ["Nuages opaques"], loc=1)
		plt.legend([p1,p2,p3,p4], ["Zone de transition","Vrais eteints","Faux eteints","Nuages"], loc=1)
		plt.title(region + ' Du ' + fichier[0][17:27] +' au ' + fichier[nbr_fichier -1][17:27] + ' Nombre de fichiers :' + str(nbr_fichier))
		plt.xlabel("Frequence d'occurence")
		plt.ylabel("Altitude en  km")
		path = region
		path += '.png'
		plt.savefig(path)
		plt.close()

def affiche_freq_separee(alt,op,cl,vr,fa,region,fichier,nbr_fichier):
	fig = plt.figure(figsize=(9, 15), dpi=100)
	p1, = plt.plot(op,alt)
	p2, = plt.plot(cl,alt)
	p3, = plt.plot(vr,alt)
	plt.legend([p1,p2,p3], ["Zone de transition","Nuages","Vrai eteint"], loc=1)
	plt.title(region + ' Du ' + fichier[0][17:27] +' au ' + fichier[nbr_fichier-1][17:27] + ' Nombre de fichiers :' + str(nbr_fichier))
	plt.xlabel("Frequence d'occurence")
	plt.ylabel("Altitude en  km")
	path = region
	path += 'occurence_zone_de_transition_nuage_separ.png'
	plt.savefig(path)
	plt.close()
	
	p3, = plt.plot(vr,alt)
	p4, = plt.plot(fa,alt)
	plt.legend([p3,p4], ["Vrais eteints","Faux eteints"], loc=1)
	plt.title(region + ' ' +str(nbr_fichier) + ' fichiers')
	plt.xlabel("Frequence d'occurence")
	plt.ylabel("Altitude en  km")
	path = region
	path += 'occurence_vrai_faux_separ.png'
	plt.savefig(path)
	plt.close()



def histo_1D_orbite(dossier,test):
	nbr_fichier,fichier = lecture_dossier(dossier)	
	nbr_fichier = int(nbr_fichier )
	lg, lt, ti, atb,alt, sol_alt, SR = extraction(dossier + fichier[0]) 
	u,v = 60000,40 # passe manuellement car certains SR ont des tailles inferieur a 56295
	print np.shape(alt)
	alt_array = np.zeros((u,v))
	alt_array[:][:] = alt
	
	if test == 1:
		for i in range(0,nbr_fichier):	
			fichier3 = dossier + fichier[i]
			print fichier3
			lg, lt, ti, atb,alt, sol_alt, SR = extraction(fichier3)
			u,v = np.shape(SR)
		 	alt_array_copy = alt_array[0:u,0:v]
			alt_max = alt_max_orbite_flag(alt,SR)
			alt_SR5 = alt_min_SR5(alt_array_copy,SR,alt_max)
			idx = ( alt_SR5 > 0 )
			temp = alt_SR5[idx]-alt_max[idx,0]
			if i == 0: total = temp
			else : total = np.concatenate([total, temp])
			print np.shape(total)
		output = '1D_distri_total'
		np.savez(output, total)

	if test == 0:
		fig = plt.figure(figsize=(14, 7), dpi=100)
		#(mu, sigma) = norm.fit(total)
		fil = '1D_distri_total.npz'
		npz = np.load(fil)
		tab = npz['arr_0']
		print tab
		mu = tab.mean()
		sigma = np.std(tab)
		plt.hist(tab, bins=np.arange(0,14,0.5))
		plt.title('Du ' + fichier[0][17:27] +' au ' + fichier[nbr_fichier-1][17:27] + ' Nombre de fichiers :' + str(nbr_fichier)+ ' \n Moyenne : ' + str(mu) + ' Ecart type :' + str(sigma))
		plt.xlabel('Taille de la zone de transition en km')
		plt.ylabel("Nombre d'occurences")
		path = 'globale_distri_taille_transition.png'
		plt.savefig(path)
		plt.close()

def histo_2D(alt, SR,lt,nbr_profil,profil_debut):
	n1 = profil_debut
	n2 = profil_debut + nbr_profil -1
	u,v = np.shape(SR)	
	alt_array = np.zeros_like(SR)
	alt_array[:][:] = alt 
	fig = plt.figure(figsize=(14, 7))
	SRbin=np.arange(-1,40)
	altbin = np.arange(-0.240,18.71+0.240,0.480)
	x = np.concatenate([SR[n1,:], SR[n2,:]], axis=0)
	y = np.concatenate([alt_array[n1,:], alt_array[n2,:]], axis=0)
	#x = SR.ravel()
	#y = alt_array.ravel()
	print 'hehe',np.shape(x),n1,n2, x

	plt.hist2d(x, y,bins=(SRbin,altbin))
	if nbr_profil == 2 :		
		plt.plot(SR[n1,:], alt_array[0,:],'-*k',linewidth=2)
		plt.plot(SR[n2,:], alt_array[0,:],'-+y',linewidth=2)
	plt.xlim(-5,40)
	plt.xlabel('SR')
	#plt.clim(0,300)
	plt.ylabel('Altitude en km')
	if nbr_profil ==1 : plt.title('Latitude :'+ str( lt[n1]))
	else : plt.title('Latitude '+ str(lt[n1])+' a '+ str(lt[n2]))
	plt.colorbar()
	plt.show()

def histo_2D_orbite(alt, SR,o,p, test):	
	alt_array = np.zeros_like(SR)
	alt_array[:][:] = alt 
	u, v = np.shape(SR)
	if test == 0:
		n1 = 0
		n2 =u * (v)	
	elif test == 1:
		n1 = o * (v)
		n2 = p * (v)

	alt_max = alt_max_orbite_flag(alt, SR)
	alt_minSR =alt_min_SR5_orbite(alt_array,SR,alt_max)

	idx = (alt_max != -10)#(alt_array <= alt_max) & (alt_array >= alt_minSR) & (alt_minSR != -1)
	SR = SR[idx]
	alt_array = alt_array[idx]
	print 'SR',np.shape(SR)
	fig = plt.figure(figsize=(14, 7))
	#SRb=np.arange(-1,40,1)
	SRb=np.concatenate([np.arange(-1,3,0.2),np.arange(3,40,1)])
	binwidth = np.max(np.shape(SRb))
	altbin = np.arange(-0.240,18.71+0.240,0.480)
	#altbin = np.linspace(-0.240,19,binwidth)
	#x = np.concatenate([SR[0,:], SR[u-1,:]], axis=0)
	#y = np.concatenate([alt_array[0,:], alt_array[u-1,:]], axis=0)
	x = SR.ravel() 
	y = alt_array.ravel()
	plt.hist2d(x[n1:n2], y[n1:n2],bins=(SRb,altbin))
	plt.xlabel('SR')
	plt.clim(0,800)
	plt.ylabel('Altitude en km')
	titre = ('Histogramme 2D de la latitude '+str(lt[o])+' a '+str(lt[p]))
	plt.title(titre)
	plt.colorbar()
	plt.show()

def histo_1D_region2(dossier,test):
	nbr_fichier,fichier = lecture_dossier(dossier)	
	nbr_fichier = int(nbr_fichier)	
	lg, lt, ti, atb,alt, sol_alt, SR = extraction(dossier + fichier[0]) 
	u,v = 60000,40 # passe manuellement car certains SR ont des tailles inferieur a 56295
	print np.shape(alt)
	alt_array = np.zeros((u,v))
	alt_array[:][:] = alt
	latbands = {'tropique' :[-20,20],'pole_nord' :[60,90],'pole_sud' : [-90,-60],'mid_sud' :[-50,-30],'mid_nord' :[30, 50 ]}
	dic_hist = dict()
	jour = np.zeros((nbr_fichier),dtype = int)
	if test == 1:
		for i in range(0,nbr_fichier):	
			fichier3 = dossier + fichier[i]
			jour[i] = int(fichier[i][25:27])
			print jour[i], fichier[i]
			lg, lt, ti, atb,alt, sol_alt, SR = extraction(fichier3)
			u,v = np.shape(SR)
		 	alt_array_copy = alt_array[0:u,0:v]
			print np.shape(alt_array_copy),np.shape(SR)
			print i, '/', nbr_fichier
			
			#idx1 = (lt < latbands[region][1]) & (lt > latbands[region][0])
			alt_max = alt_max_orbite_flag(alt,SR)
			alt_SR5 = alt_min_SR5(alt_array_copy,SR,alt_max)
			idx = ( alt_SR5 > 0 ) & (lt < 20) & (lt > -20)
			temp = alt_SR5[idx]-alt_max[idx,0]
			if i == 0: dic_hist['tropique'] = temp
			else : dic_hist['tropique'] = np.concatenate([dic_hist['tropique'], temp])

			idx = ( alt_SR5 > 0 ) & (lt < 90) & (lt > 60)
			temp = alt_SR5[idx]-alt_max[idx,0]
			if i == 0: dic_hist['pole_nord'] = temp
			else : dic_hist['pole_nord'] = np.concatenate([dic_hist['pole_nord'], temp])

			idx = ( alt_SR5 > 0 ) & (lt < -60) & (lt >-90)
			temp = alt_SR5[idx]-alt_max[idx,0]
			if i == 0: dic_hist['pole_sud'] = temp
			else : dic_hist['pole_sud'] = np.concatenate([dic_hist['pole_sud'], temp])

			idx = ( alt_SR5 > 0 ) & (lt < -30) & (lt > -50)
			temp = alt_SR5[idx]-alt_max[idx,0]
			if i == 0: dic_hist['mid_sud'] = temp
			else : dic_hist['mid_sud']  = np.concatenate([dic_hist['mid_sud'] , temp])		

			idx = ( alt_SR5 > 0 ) & (lt < 50) & (lt > 30)
			temp = alt_SR5[idx]-alt_max[idx,0]
			if i == 0: dic_hist['mid_nord'] = temp
			else : dic_hist['mid_nord']  = np.concatenate([dic_hist['mid_nord'] , temp])
		for	region in latbands.keys():
			output = region +'1D_distri'
			np.savez(output, dic_hist[region])

	if test == 0:
		for region in latbands.keys():
			fil = region +'1D_distri.npz'
			npz = np.load(fil)
			tab = npz['arr_0']
			mu = tab.mean()
			sigma = np.std(tab)
			plt.hist(tab, bins=np.arange(0,14,0.5))
			plt.rcParams["axes.titlesize"] = 8
			plt.title(region+'  Du ' + fichier[0][17:25]+str(np.min(jour)) +' au ' + fichier[nbr_fichier-1][17:25] +str(np.max(jour))+ ' Nombre de fichiers :' + str(nbr_fichier)+ '\n Moyenne : '+ str(mu) +' Ecart type : '+ str(sigma))
			plt.xlabel('Taille de la zone de transition en km')
			plt.ylabel("Nombre d'occurences")
			path = region
			path += '_histo_1D_moins.png'
			plt.savefig(path)
			plt.close()
			#plt.show()

def histo_2D_region(SR,alt,dossier):
	nbr_fichier,fichier = lecture_dossier(dossier)	
	nbr_fichier = int(nbr_fichier / 1.)
	latbands = {'tropique' :[-20,20],'pole_nord' :[60,90],'pole_sud' : [-90,-60],'mid_sud' :[-50,-30],'mid_nord' :[30, 50 ]}
	dic_hist = dict()
	jour = np.zeros((nbr_fichier),dtype = int)
	print nbr_fichier
	u,v = 100000*(nbr_fichier)/5,40 # passe manuellement car certains SR ont des tailles inferieur a 56295
	print np.shape(alt)
	alt_array = np.zeros((u,v))
	alt_array[:][:] = alt
	alt_array = alt_array.ravel()
	SRb=[0,0.01,1.2,3,5,7,10,15,20,25,30,40,50,60,80]
	altbin = np.arange(-0.240,18.71+0.240,0.480)

	for i in range(0,nbr_fichier):	
		fichier3 = dossier + fichier[i]
		jour[i] = int(fichier[i][25:27])
		print jour[i], fichier[i]
		lg, lt, ti, atb,alt, sol_alt, SR = extraction(fichier3)
		print i, '/', nbr_fichier
		
		#idx1 = (lt < latbands[region][1]) & (lt > latbands[region][0])
	
		idx = (lt < 20) & (lt > -20)
		SR_copy = SR[idx,:].ravel()
		o = np.max(np.shape(SR_copy))
		alt_array_copy = alt_array[0:o]	
		temp, xx , yy = np.histogram2d(SR_copy, alt_array_copy, bins = (SRb, altbin))
		if i == 0: dic_hist['tropique'] = temp
		else : dic_hist['tropique'] = dic_hist['tropique'] + temp

		idx = (lt < 90) & (lt > 60)
		SR_copy = SR[idx,:].ravel()
		o = np.max(np.shape(SR_copy))
		alt_array_copy = alt_array[0:o]	
		temp, xx , yy = np.histogram2d(SR_copy, alt_array_copy, bins = (SRb, altbin)) 
		if i == 0: dic_hist['pole_nord'] = temp
		else : dic_hist['pole_nord'] = dic_hist['pole_nord'] + temp

		idx = (lt < -60) & (lt > -90)
		SR_copy = SR[idx,:].ravel()
		o = np.max(np.shape(SR_copy))
		alt_array_copy = alt_array[0:o]	
		temp, xx , yy = np.histogram2d(SR_copy, alt_array_copy, bins = (SRb, altbin))
		if i == 0: dic_hist['pole_sud'] = temp
		else : dic_hist['pole_sud'] = dic_hist['pole_sud'] + temp

		idx = (lt < -30) & (lt > -50)
		SR_copy = SR[idx,:].ravel()
		o = np.max(np.shape(SR_copy))
		alt_array_copy = alt_array[0:o]	
		temp, xx , yy = np.histogram2d(SR_copy, alt_array_copy, bins = (SRb, altbin)) 
		if i == 0: dic_hist['mid_sud'] = temp
		else : dic_hist['mid_sud'] = dic_hist['mid_sud'] + temp	

		idx = (lt < 50) & (lt > 30)
		SR_copy = SR[idx,:].ravel()
		o = np.max(np.shape(SR_copy))
		alt_array_copy = alt_array[0:o]	
		temp, xx , yy = np.histogram2d(SR_copy, alt_array_copy, bins = (SRb, altbin))
		if i == 0: dic_hist['mid_nord'] = temp
		else : dic_hist['mid_nord'] = dic_hist['mid_nord'] + temp


		SR_copy = SR.ravel()
		o = np.max(np.shape(SR_copy))
		alt_array_copy = alt_array[0:o]	
		temp, xx , yy = np.histogram2d(SR_copy, alt_array_copy, bins = (SRb, altbin))
		if i == 0: dic_hist['total'] = temp
		else : dic_hist['total'] = dic_hist['total'] + temp

	for region, tab in dic_hist.items():
		x = np.arange(0,np.shape(SRb)[0])
		plt.pcolormesh(x, altbin, dic_hist[region].T / np.sum(dic_hist[region]) )
		plt.gca().set_xticks(x)
		plt.gca().set_xticklabels(SRb)

		plt.rcParams["axes.titlesize"] = 10
		plt.title(region+'  Du ' + fichier[0][17:25]+str(np.min(jour)) +' au ' + fichier[0][17:25] +str(np.max(jour))+ ' Nombre de fichiers :' + str(nbr_fichier))
		plt.clim(0,0.001)
		#plt.xlim(0,80)
		plt.ylim(0,20)

		plt.xlabel('SR')
		plt.ylabel("Altitude en km")
		plt.colorbar()
		path = region
		path += '_histo_2D.png'
		plt.savefig(path)
		plt.clf()

def hist_cloud_ext(fichier,date):
	lg, lt, ti, atb,alt, sol_alt, SR = extraction(fichier)
	SRb=[0,0.01,1.2,3,5,7,10,15,20,25,30,40,50,60,80]
	altbin = np.arange(-0.240,18.71+0.240,0.480)
	u,v =  np.shape(SR)
	alt_array = np.zeros((u,v))
	alt_array [:][:] = alt
	alt_max = alt_max_orbite_flag(alt, SR)
	#idx = ( alt_max != -10)

	alt_max = alt_max_orbite_flag(alt, SR)
	alt_minSR =alt_min_SR5_orbite(alt_array,SR,alt_max)

	idx = (alt_max != -10) & ((alt_array <= alt_max) | (alt_array >= alt_minSR)) & (alt_minSR != -1)
	
	SR = SR[idx].ravel()
	alt_array= alt_array[idx].ravel()
	fig = plt.figure(figsize=(14, 7), dpi=100)
	#plt.hist2d(SR,alt_array,bins= (SRb,altbin))
	histSR, x,y = np.histogram2d(SR ,alt_array,bins = (SRb, altbin))
	x = np.arange(0,np.shape(SRb)[0])
	plt.pcolormesh(x, altbin, histSR.T / np.sum(histSR))
	plt.gca().set_xticks(x)
	plt.gca().set_xticklabels(SRb)
	plt.xlabel('SR')
	plt.ylabel('Altitude en km')
	plt.ylim(0,19)
	plt.title('Histogramme 2D sur une orbite des nuages responsables des extinctions')
	plt.clim(0,0.001)
	plt.colorbar()
	path = './lot/'+ date + '.png'
	plt.show()
	plt.savefig(path)
	plt.clf
	plt.close()

def hist_could_ext_dossier(dossier):
	nbr_fichier,fichier = lecture_dossier(dossier)	
	nbr_fichier = int( nbr_fichier )	
	for i in range(0,nbr_fichier):
		print '== fichier', fichier[i]
		print i ,'/', nbr_fichier
		calfile = dossier + fichier[i]
		date = fichier[i][17:27]
		hist_cloud_ext(calfile,date)



def histo_2D_region_eteint(SR,alt,dossier,test):
	nbr_fichier,fichier = lecture_dossier(dossier)	
	nbr_fichier = int(nbr_fichier)
	latbands = {'tropique' :[-20,20],'pole_nord' :[60,90],'pole_sud' : [-90,-60],'mid_sud' :[-50,-30],'mid_nord' :[30, 50 ]}
	dic_hist = dict()
	jour = np.zeros((nbr_fichier),dtype = int)
	print nbr_fichier
	u,v = 100000*(nbr_fichier)/5,40 # passe manuellement car certains SR ont des tailles inferieur a 56295
	print np.shape(alt)
	alt_array = np.zeros((u,v))
	alt_array[:][:] = alt
	alt_array = alt_array.ravel()
	SRb=[0,0.01,1.2,3,5,7,10,15,20,25,30,40,50,60,80]
	altbin = np.arange(-0.240,18.71+0.240,0.480)

	if test == 1: 
		for i in range(0,nbr_fichier):	
			fichier3 = dossier + fichier[i]
			jour[i] = int(fichier[i][25:27])
			print jour[i], fichier[i]
			lg, lt, ti, atb,alt, sol_alt, SR = extraction(fichier3)
			alt_max = alt_max_orbite(alt, SR)
			print i, '/', nbr_fichier
			
			#idx1 = (lt < latbands[region][1]) & (lt > latbands[region][0])
		
			idx = (lt < 20) & (lt > -20) & (alt_max != -10  )
			SR_copy = SR[idx,:].ravel()
			o = np.max(np.shape(SR_copy))
			alt_array_copy = alt_array[0:o]	
			temp, xx , yy = np.histogram2d(SR_copy, alt_array_copy, bins = (SRb, altbin))
			if i == 0: dic_hist['tropique'] = temp
			else : dic_hist['tropique'] = dic_hist['tropique'] + temp

			idx = (lt < 90) & (lt > 60) & (alt_max != -10  )
			SR_copy = SR[idx,:].ravel()
			o = np.max(np.shape(SR_copy))
			alt_array_copy = alt_array[0:o]	
			temp, xx , yy = np.histogram2d(SR_copy, alt_array_copy, bins = (SRb, altbin)) 
			if i == 0: dic_hist['pole_nord'] = temp
			else : dic_hist['pole_nord'] = dic_hist['pole_nord'] + temp

			idx = (lt < -60) & (lt > -90) & (alt_max != -10  )
			SR_copy = SR[idx,:].ravel()
			o = np.max(np.shape(SR_copy))
			alt_array_copy = alt_array[0:o]	
			temp, xx , yy = np.histogram2d(SR_copy, alt_array_copy, bins = (SRb, altbin))
			if i == 0: dic_hist['pole_sud'] = temp
			else : dic_hist['pole_sud'] = dic_hist['pole_sud'] + temp

			idx = (lt < -30) & (lt > -50) & (alt_max != -10  )
			SR_copy = SR[idx,:].ravel()
			o = np.max(np.shape(SR_copy))
			alt_array_copy = alt_array[0:o]	
			temp, xx , yy = np.histogram2d(SR_copy, alt_array_copy, bins = (SRb, altbin)) 
			if i == 0: dic_hist['mid_sud'] = temp
			else : dic_hist['mid_sud'] = dic_hist['mid_sud'] + temp	

			idx = (lt < 50) & (lt > 30) & (alt_max != -10  )
			SR_copy = SR[idx,:].ravel()
			o = np.max(np.shape(SR_copy))
			alt_array_copy = alt_array[0:o]	
			temp, xx , yy = np.histogram2d(SR_copy, alt_array_copy, bins = (SRb, altbin))
			if i == 0: dic_hist['mid_nord'] = temp
			else : dic_hist['mid_nord'] = dic_hist['mid_nord'] + temp
			
		for region in latbands.keys():
			output = region+'_2D_extok'
			np.savez(output,dic_hist[region],alt)


	if test == 0:

		for region in latbands.keys():
			fil = region+'_2D_extok.npz'
			npz = np.load(fil)
			tab = npz['arr_0']

			x = np.arange(0,np.shape(SRb)[0])
			plt.pcolormesh(x, altbin, tab.T / np.sum(tab) )
			plt.gca().set_xticks(x)
			plt.gca().set_xticklabels(SRb)

			plt.rcParams["axes.titlesize"] = 10
			plt.title(region + ' Eteints' +'  Du ' + fichier[0][17:25]+str(np.min(jour)) +' au ' + fichier[0][17:25] +str(np.max(jour))+ ' Nombre de fichiers :' + str(nbr_fichier))
			plt.clim(0,0.001)
			#plt.xlim(0,80)
			plt.ylim(0,20)

			plt.xlabel('SR')
			plt.ylabel("Altitude en km")
			plt.colorbar()
			path = region
			path += '_histo_2D_eteint.png'
			
			plt.savefig(path)
			plt.clf()

def histo_2D_warmpool_eteint(SR,alt,dossier,test):
	nbr_fichier,fichier = lecture_dossier(dossier)	
	nbr_fichier = int(nbr_fichier)
	jour = np.zeros((nbr_fichier),dtype = int)
	print nbr_fichier
	u,v = 100000*(nbr_fichier)/5,40 # passe manuellement car certains SR ont des tailles inferieur a 56295
	print np.shape(alt)
	alt_array = np.zeros((u,v))
	alt_array[:][:] = alt
	alt_array = alt_array.ravel()
	SRb=[0,0.01,1.2,3,5,7,10,15,20,25,30,40,50,60,80]
	altbin = np.arange(-0.240,18.71+0.240,0.480)
	if test == 1:
		for i in range(0,nbr_fichier):	
			fichier3 = dossier + fichier[i]
			jour[i] = int(fichier[i][25:27])
			print jour[i], fichier[i]
			lg, lt, ti, atb,alt, sol_alt, SR = extraction(fichier3)
			alt_max = alt_max_orbite(alt, SR)
			print i, '/', nbr_fichier
		
		
			#idx = (lt < 15) & (lt > -15) & (alt_max != -10  ) & ( lg<180 ) & (lg >70) #warmpool
			
			idx = (lt < -25) & (lt > -45) & (alt_max != -10  ) & ( lg<-85 ) & (lg >-95) #shallowcu
			
			#idx =  (lt < 0) & (lt > -15) & (alt_max != -10  ) & ( lg<-75 ) & (lg >-90) #stratocu
			SR_copy = SR[idx,:].ravel()
			o = np.max(np.shape(SR_copy))
			alt_array_copy = alt_array[0:o]	
			temp, xx , yy = np.histogram2d(SR_copy, alt_array_copy, bins = (SRb, altbin))
			if i == 0: total = temp
			else : total = total + temp

		output = 'shlcu'
		np.savez(output,total)
	if test == 0:
		fil = 'shlcu.npz'
		npz = np.load(fil)
		total = npz['arr_0']
		x = np.arange(0,np.shape(SRb)[0])
		plt.pcolormesh(x, altbin, total.T / np.sum(total) )
		plt.gca().set_xticks(x)
		plt.gca().set_xticklabels(SRb)

		plt.rcParams["axes.titlesize"] = 10
		plt.title('Shallowcu ' +'  Du ' + fichier[0][17:25]+str(np.min(jour)) +' au ' + fichier[0][17:25] +str(np.max(jour))+ ' Nombre de fichiers :' + str(nbr_fichier))
		plt.clim(0,0.001)
		#plt.xlim(0,80)
		plt.ylim(0,20)

		plt.xlabel('SR')
		plt.ylabel("Altitude en km")
		plt.colorbar()
		#path = '_warmpool.png'
		#path = '_stratocu.png'
		path = '_shallow_vrai.png'
		plt.savefig(path)
		plt.clf()

def distri_alt_max(dossier,test) :
	
	dic_distri = dict()
	latbands = {'tropique' :[-20,20],'pole_nord' :[60,90],'pole_sud' : [-90,-60],'mid_sud' :[-50,-30],'mid_nord' :[30, 50 ]}
	nbr_fichier,fichier = lecture_dossier(dossier)	
	nbr_fichier = int(nbr_fichier )
	if test == 1 :
		for i in range(0,nbr_fichier):
			print i, '/',nbr_fichier
			lg, lt, ti, atb,alt, sol_alt, SR = extraction(dossier + fichier[i])
			alt_max = alt_max_orbite(alt, SR)
			for region in latbands.keys():
				idx = (lt < latbands[region][1]) & (lt > latbands[region][0])
				temp = alt_max[idx]
				if i == 0 : dic_distri[region] = temp
				else : dic_distri[region] = np.concatenate([dic_distri[region],temp])
		for	region in latbands.keys():
			output = region + 'dst_max'
			np.savez(output, dic_distri[region],alt)

	if test == 0 :

		for region in latbands.keys():
			fil = region + 'dst_max'+'.npz'
			npz = np.load(fil)
			tab = npz['arr_0']
			alt= npz['arr_1']
			tab = np.ma.masked_where(tab == -10, tab)
			bin = np.arange(0,14,0.5)
			#plt.hist(hist[0],hist[1],orientation = 'horizontal' )
			idx = tab >=0
			me = tab[idx]
			plt.hist(tab, bins=bin,orientation = 'horizontal')
			plt.rcParams["axes.titlesize"] = 8
			mu = me.mean()
			sigma = np.std(me)
			plt.title(region+' Nombre de fichiers : '+str(nbr_fichier) + '\n Moyenne :' + str(mu) +' km Ecart type : ' + str(sigma) +' km')
			plt.xlabel('Nombre d occurences')
			plt.ylabel("Epaisseur de la zone de transition (km)")
			plt.ylim(0,19)
			path = region
			path += 'distri_alt_max.png'
			
			plt.savefig(path)
			plt.close()


def hist_lonlat(dossier,test):		
	nbr_fichier,fichier = lecture_dossier(dossier)	
	nbr_fichier = int(nbr_fichier)

	if test ==1 :

		for i in range(0,nbr_fichier):
			lg, lt, ti, atb,alt, sol_alt, SR = extraction(dossier + fichier[i])
			print np.shape(lt),np.shape(lg)
			alt_max = alt_max_orbite(alt,SR)
			idx = (alt_max != -10 )
			alt_max = alt_max[idx]
			lg = lg[idx]	
			lt = lt[idx]
			lonb = np.arange(-90,90,2)
			latb = np.arange(-180,180,2)
			if i == 0:
				latotal = lt
				lgtotal = lg
			else:
				latotal = np.concatenate([latotal,lt])
				lgtotal = np.concatenate([lgtotal,lg])

		output = 'lonlat_total'
		np.savez(output,latotal,lgtotal)

	if test ==0:	
		fig = plt.figure(figsize=(14, 7), dpi=100)
		m = Basemap()
		fil = 'lonlat_total.npz'
		npz = np.load(fil)
		latotal = npz['arr_0']
		lgtotal = npz['arr_1']
		binlt = np.arange(-90,90,4)
		binlg = np.arange(-180,180,4)
		hist,xx,yy = np.histogram2d( lgtotal, latotal, bins = (binlg,binlt) )
		#x = np.arange(0,np.shape(SRb)[0])
		m.drawcoastlines(color ='white')
		plt.pcolormesh(binlg, binlt, 100.* hist.T / np.sum(hist),alpha = 1 )
		plt.rcParams["axes.titlesize"] = 8
		plt.title('  Distribution de alt max - Nombre de fichiers : '+str(nbr_fichier))
		plt.colorbar()
		#m.bluemarble()
		path = 'lon_lat2.png'
		plt.xlim(-180,180)
		plt.ylim(-90,90)
		plt.savefig(path)
		#plt.show()
		plt.close()	

	
if __name__=='__main__':

	fichier3="/bdd/CFMIP/CFMIP_OBS_LOCAL/GOCCP/instant_SR_CR_DR/grid_L40/2013/201308/night/instant_SR_CR_DR_2013-08-24T23-19-30ZN_night_CFMIP2_2.70.nc"#38418 mid
	#fichier3="/bdd/CFMIP/CFMIP_OBS_LOCAL/GOCCP/instant_SR_CR_DR/grid_L40/2013/201308/night/instant_SR_CR_DR_2013-08-23T19-18-36ZN_night_CFMIP2_2.70.nc" #20632 trop
	listevar(fichier3)
	n1=0
	n1 , n2 = 0  , 20632
	lg, lt, ti, atb,alt, sol_alt, SR = extraction(fichier3)
 	#alt_max(SR, alt, n1,2*n2)
	#profil_SR(n1,2*n2,lg,lt,ti,SR,alt)
	#affiche(n1,n2,lg, lt, ti, SR,alt,sol_alt)
	'''
	n1 = 0
	n2, lt, lg, ti, atb, sol_alt = coord_pole(lg, lt, ti, atb,sol_alt)
	alt_max(SR, alt, n1,n2)
	#profil(n1,n2,lg,lt,ti,atb,alt)	
	affiche(n1,n2,lg, lt, ti, atb,alt,sol_alt)
	lonlat(n1, n2,lg, lt, ti, atb)
	'''

	#lonlat(n1, n2,lg, lt, ti, atb)
	#affiche(n1,n2,lg, lt, ti, SR,alt,sol_alt)
	#alt_max_orbite(alt, sol_alt, SR,lt)
	
	# Tableau de flag
	
	#affiche_flag(lg, lt, flag,alt,sol_alt)
	#freq_occur(flag,alt)
	#n2,x = np.shape(SR)
	dossier = '/bdd/CFMIP/CFMIP_OBS_LOCAL/GOCCP/instant_SR_CR_DR/grid_L40/2013/201308/night/'
	#freq_occur_region(dossier)

	#lonlat(n1, n2,lg, lt, ti, atb)
	afficheSR(lg, lt, ti, SR,alt,sol_alt)
	n = 38418
	#profil_SR(n,lg,lt,ti,SR,alt,fichier3)
	#n = 41750
	#profil(n,lg,lt,ti,atb,alt,fichier3)
	#freq_occur_region(dossier,0)
	#alti = alt_min_SR5(alt,SR)
	#flag = flag_orbite(alt,SR)
	#affiche_flag(lg, lt, flag,alt,sol_alt)
	#n1,n2 =  38352 , 42466
	#histo_2D_orbite(alt,SR,n1,n2,0)#test : 1 = on definit les bornes 0 = sur l'orbite
	#alt_min_SR5(alt,SR)
	#histo_1D_region2(dossier,0)
	#histo_1D_orbite(dossier,0)
	#histo_1D(SR,alt)
	#nbr = 2
	#debut = 38418
	#histo_2D(alt,SR,lt,nbr,debut)
	#histo_2D_region_eteint(SR,alt,dossier,0)	
	#hist_cloud_ext(fichier3,'test_normalisation')
	#hist_could_ext_dossier(dossier)
	#distri_alt_max(dossier,0)
	histo_2D_warmpool_eteint(SR,alt,dossier,1)
	#hist_lonlat(dossier,0)
	#draw()
	#listevar(fichier3)