1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185
| intervalles=list()
intervalles.append(100)
intervalles.append(50)
intervalles.append(25)
intervalles.append(10)
pourcentage_liste=list()
pourcentage_liste.append(10)
pourcentage_liste.append(5)
pourcentage_liste.append(2.5)
pourcentage_listee.append(1)
for i in intervalles:
for k in pourcentage_liste:
(pts_inf, pts_sup), (pts_inf_mean, pts_sup_mean) = courbes_phi(fvc, temp_surf, nb_interval=i, pourcentage=k)
# Liste de couleurs pour l'affichage des points des intervalles
nb_interval=i
pourcentage=k
color_list = plt.cm.Set1(np.linspace(0, 1, nb_interval))
fig, ax = plt.subplots(figsize=(14, 10))
fig.suptitle("Estimation du coefficient de Priestley-Taylor \n par double interpolation", fontsize=20)
x_reg = np.array([0, 1])
# Points moyens
print("*** Régression avec les points moyens ***")
reg_h_moy = linregress(pts_inf_mean[:, 0], pts_inf_mean[:, 1])
print("Bord humide : R2 = %f - pente = %f - intersect = %f" %
(reg_h_moy[2] ** 2, reg_h_moy[0], reg_h_moy[1]))
y_inf = reg_h_moy[0]*x_reg + reg_h_moy[1]
# Affichage des points moyens du bord humide
ax.scatter(pts_inf_mean[:, 0], pts_inf_mean[:, 1], zorder=10, c="k")
# Affichage de la régression des points moyens du bord humide
ax.plot(x_reg, y_inf,'k-')
reg_s_moy = linregress(pts_sup_mean[:, 0], pts_sup_mean[:, 1])
print("Bord sec : R2 = %f - pente = %f - intersect = %f" %
(reg_s_moy[2] ** 2, reg_s_moy[0], reg_s_moy[1]))
y_sup = reg_s_moy[0]*x_reg + reg_s_moy[1]
# Affichage des points moyens du bord sec
pts = ax.scatter(pts_sup_mean[:, 0], pts_sup_mean[:, 1], c="k", zorder=10, label="Points moyens")
# Affichage de la régression des points moyens du bord sec
reg_moy, = ax.plot(x_reg, y_sup,'k-', label="Régression des points moyens")
# Affichage doite du bord humide
temp_hum = pts_inf_mean.copy()
arg = np.argsort(temp_hum[:, 1])
y_min_hum = temp_hum[:, 1][arg][1]
plot_reg_humide, = ax.plot(x_reg, [y_min_hum, y_min_hum],'r-', lw=2, label="Droite du bord humide")
# Tous les points
print("*** Régression avec l'ensemble des points ***")
reg_h_tot = linregress(pts_inf[:, 0], pts_inf[:, 1])
print("Bord humide : R2 = %f - pente = %f - intersect = %f" %
(reg_h_tot[2] ** 2, reg_h_tot[0], reg_h_tot[1]))
slope_sec = reg_h_tot[0]
y_inf = reg_h_tot[0]*x_reg + reg_h_tot[1]
# Affichage de tous les points du bord humide
ax.scatter(pts_inf[:, 0], pts_inf[:, 1], color=color_list[pts_inf[:, 2].astype(int)],
alpha=0.2, edgecolor='none')
# Affichage de la régression de tous les points du bord humide
ax.plot(x_reg, y_inf,'b-')
reg_s_tot = linregress(pts_sup[:, 0], pts_sup[:, 1])
print("Bord sec : R2 = %f - pente = %f - intersect = %f" %
(reg_s_tot[2] ** 2, reg_s_tot[0], reg_s_tot[1]))
y_sup = reg_s_tot[0]*x_reg + reg_s_tot[1]
# Affichage de tous les points du bord sec
ax.scatter(pts_sup[:, 0], pts_sup[:, 1], color=color_list[pts_sup[:, 2].astype(int)],
alpha=0.2, edgecolor="none")
# Affichage de la régression de tous les points du bord sec
plot_reg_tot, = ax.plot(x_reg, y_sup,'b-', label="Régression de l'ensemble des points")
ax.set_title("\nNombre d'intervalle = %d - Poucentage = %.1f%%" % (nb_interval, pourcentage), fontsize=16)
ax.set_xlabel("FVC", fontsize=16)
ax.set_ylabel("Ts", fontsize=16)
plt.legend(handles=[pts, reg_moy, plot_reg_tot, plot_reg_humide])
# Sauvegarde de la figure au format png
plt.savefig(os.path.join(path_indice, "phi_%d_%.1f.png" % (nb_interval, pourcentage)), dpi=300)
plt.show()
#ESTIMER PHI (COEFF DE PRIESTLEY-TAYLOR)
t_min_hum = y_min_hum
slope_sec = -18.813213
intercept_sec = 3101.330539
temp_max = slope_sec * fvc + intercept_sec
stats_temp_max=get_stats(temp_max)
#hist(temp_max)
#create_tif(temp_max, "temp_max_%d_%.1f.tif" % (nb_interval, pourcentage))
phi = 1.26 * (temp_max - temp_surf ) / (temp_max - t_min_hum)
stats_phi=get_stats(phi)
#hist(phi)
#create_tif(phi, "phi_%d_%.1f_ts.tif" % (nb_interval, pourcentage))
print("Nombre de valeurs négative : %d" % (phi < 0).sum())
print("Nombre de valeurs > 1,26 : %d" % (phi > 1.26).sum())
val_negative=((phi < 0).sum())
val_sup_1_26=((phi > 1.26).sum())
#CALCULER EF
tmoy = temp_surf
delta = (get_e_sat(tmoy) / tmoy) #* ((6790.4985 / tmoy) - 5.02808)
stats_delta=get_stats(delta)
#hist(delta)
#create_tif(delta, "delta_%d_%.1f.tif" % (nb_interval, pourcentage))
# gamma : constante psychométrique
gamma = 66
EF_WDI = delta * phi / (delta + gamma)
stats_EF=get_stats(EF_WDI)
#hist(EF_WDI)
#create_tif(EF_WDI, "EF_%d_%.1f_ts.tif" % (nb_interval, pourcentage))
#STOCKAGE DES STATS
res_1_tour=list()
#méthode de calcul de EF
res_1_tour.append("WDI")
#nombre d intervalle
res_1_tour.append(i)
#pourcentage de points
res_1_tour.append(k)
#pour les points moyens
#R2
res_1_tour.append(reg_h_moy[2] ** 2)
#pente
res_1_tour.append(reg_h_moy[0])
#intercept
res_1_tour.append(reg_h_moy[1])
#R2
res_1_tour.append(reg_s_moy[2] ** 2)
#pente
res_1_tour.append(reg_s_moy[0])
#intercept
res_1_tour.append(reg_s_moy[1])
#pour l'ensemble des points
#R2
res_1_tour.append(reg_h_tot[2] ** 2)
#pente
res_1_tour.append(reg_h_tot[0])
#intercept
res_1_tour.append(reg_h_tot[1])
#R2
res_1_tour.append(reg_s_tot[2] ** 2)
#pente
res_1_tour.append(reg_s_tot[0])
#intercept
res_1_tour.append(reg_s_tot[1])
#pour la temperature max
#valeur min
res_1_tour.append(stats_temp_max[0])
#valeur moy
res_1_tour.append(stats_temp_max[2])
#valeur max
res_1_tour.append(stats_temp_max[1])
#pour le coefficient de priest-taylor(phi)
#nombre de valeur negative
res_1_tour.append(val_negative)
#nombre de valeur superieur a 1.26
res_1_tour.append(val_sup_1_26)
#valeur min
res_1_tour.append(stats_phi[0])
#valeur moy
res_1_tour.append(stats_phi[2])
#valeur max
res_1_tour.append(stats_phi[2])
#pour la pente de la droite de saturation de l air en vapeur d eau (delta)
#valeur min
res_1_tour.append(stats_delta[0])
#valeur moy
res_1_tour.append(stats_delta[2])
#valeur max
res_1_tour.append(stats_delta[1])
#pour EF
#valeur min
res_1_tour.append(stats_EF[0])
#valeur moy
res_1_tour.append(stats_EF[2])
#valeur max
res_1_tour.append(stats_EF[1]) |
Partager