Changeset 51

Timestamp:

08/08/14 18:09:30 (10 years ago)

Author:

lahlod

Message:

new scripts

Location:

trunk/src/scripts_Laura/ARCTIC/Travail_CEN

Files:

• daily_emis_AMSUA_AMSUB_89.py (added)
• daily_evolution_anomaly_grad_ratio.py (added)
• daily_gradient_ratio.py (added)
• daily_param_emissivity_all_arctic.py (added)
• map_monthly_mean_emis_parameters_whole_arctic.py (added)
• map_statistic_parameters.py (modified) (1 diff)
• test.py (modified) (3 diffs)

trunk/src/scripts_Laura/ARCTIC/Travail_CEN/map_statistic_parameters.py

@@ -211 +211 @@
 # time evolution (monthly) in a given zone #
 ############################################
-# read monthly std of emis 
-ifr = 3
-fichier_temporal_std = Dataset('/net/argos/data/parvati/lahlod/ARCTIC/monthly_GLACE/gridded_data/cartesian_grid/res_40/cartesian_grid_monthly_mean-std_data_lamb_spec_near_nadir_AMSUB' + str(frequ[ifr]) + '_' + month[imo] + '2009.nc', 'r', format='NETCDF3_CLASSIC')
-
-
-
 #zone 1 : yi = 400. // yf = 720. // xi = -1080. // xf = -920. (Canadian Archipelago)
 #zone 2 : yi = 920. // yf = 1280. // xi = -720. // xf = -560. (Chukchi Sea, Bering Detroit)

trunk/src/scripts_Laura/ARCTIC/Travail_CEN/test.py

@@ -13 +13 @@
 import ffgrid2
 import map_ffgrid
+from matplotlib import colors
 from matplotlib.font_manager import FontProperties
 import map_cartesian_grid
 
 
-
-
-
-
+###############################
+# time period characteristics #
+###############################
 MONTH = np.array(['01', '02', '03', '04', '05', '06', '07', '08', '09', '10', '11', '12'])
 month = np.array(['JANUARY', 'FEBRUARY', 'MARCH', 'APRIL', 'MAY', 'JUNE', 'JULY', 'AUGUST', 'SEPTEMBER', 'OCTOBER', 'NOVEMBER', 'DECEMBER'])
@@ -27 +27 @@
 
 
-
-
-
-frequ = 23
-
-
-#########################
-# emissivity thresholds #
-#########################
-fichier = open('/net/argos/data/parvati/lahlod/ARCTIC/AMSUA_ice_class/AMSUA'+ str(frequ) + '_data_classification_parameters_ice_no-ice_2009.dat', 'r')
-numlines = 0
-for line in fichier: numlines += 1
-
-fichier.close()
-
-fichier=open('/net/argos/data/parvati/lahlod/ARCTIC/AMSUA_ice_class/AMSUA' + str(frequ) + '_data_classification_parameters_ice_no-ice_2009.dat','r')
-nbtotal = numlines-1
-iligne = 0
-mois = np.zeros([nbtotal],object)
-emis_lim_spec = np.zeros([nbtotal],float)
-emis_lim_lamb = np.zeros([nbtotal],float) 
-while (iligne < nbtotal) :
-    line=fichier.readline()
-    liste = line.split()
-    mois[iligne] = str(liste[0])
-    emis_lim_spec[iligne] = float(liste[7])
-    emis_lim_lamb[iligne] = float(liste[8])
-    iligne=iligne+1
-
-fichier.close()
-vec = np.arange(0, nbtotal + 1, 50)
-lim_spec = np.zeros([M], float)
-lim_lamb = np.zeros([M], float)
-for imo in range (0, M):
-    lim_spec[imo] = emis_lim_spec[vec[imo]]
-    lim_lamb[imo] = emis_lim_lamb[vec[imo]]
-
-
-
+########################
+# grid characteristics #
+########################
 x0 = -3000. # min limit of grid
 x1 = 2500. # max limit of grid
@@ -79 +43 @@
 
 
-emis_spec_moy = np.zeros([ny, nx, M], float)
-emis_lamb_moy = np.zeros([ny, nx, M], float)
-emis_diff = np.zeros([ny, nx, M], float)
-emis_ratio = np.zeros([ny, nx, M], float)
+
+
+grad_rat = np.zeros([M, 31, ny, nx], float)
+spec23 = np.zeros([M, 31, ny, nx], float)
+spec89 = np.zeros([M, 31, ny, nx], float)
+ratio89 = np.zeros([M, 31, ny, nx], float)
+for imo in range (0, M): 
+    # daily read gradient ratio file
+    print 'read daily gradient ratio for month ' + month[imo]
+    fichier_grad_ratio = Dataset('/net/argos/data/parvati/lahlod/ARCTIC/monthly_GLACE/gridded_data/cartesian_grid/res_40/cartesian_grid_daily_grad_ratio_spec_23-89_near_nadir_AMSUA_' + month[imo] + '2009.nc', 'r', format = 'NETCDF3_CLASSIC')
+    gr = fichier_grad_ratio.variables['grad_ratio'][:]
+    fichier_grad_ratio.close()
+    # read daily emis anomaly file for 23GHz
+    print 'read daily emis anomaly 23GHz for month ' + month[imo]
+    fichier_anom23 = Dataset('/net/argos/data/parvati/lahlod/ARCTIC/monthly_GLACE/gridded_data/cartesian_grid/res_40/cartesian_grid_daily_data_lamb_spec_ratio_anomaly_near_nadir_AMSUA23_' + month[imo] + '2009.nc', 'r', format = 'NETCDF3_CLASSIC')
+    sa_23 = fichier_anom23.variables['spec_anomaly'][:]
+    fichier_anom23.close()
+    # read daily emis anomaly file for 89GHz
+    print 'read daily emis anomaly 89GHz for month ' + month[imo]
+    fichier_anom89 = Dataset('/net/argos/data/parvati/lahlod/ARCTIC/monthly_GLACE/gridded_data/cartesian_grid/res_40/cartesian_grid_daily_data_lamb_spec_ratio_anomaly_near_nadir_AMSUA89_' + month[imo] + '2009.nc', 'r', format = 'NETCDF3_CLASSIC')
+    sa_89 = fichier_anom89.variables['spec_anomaly'][:]
+    ra_89 = fichier_anom89.variables['ratio_anomaly'][:]
+    fichier_anom89.close()
+    for ijr in range (0, month_day[imo]):
+        grad_rat[imo, ijr, :, :] = gr[ijr, :, :]
+        spec23[imo, ijr, :, :] = sa_23[ijr, :, :]
+        spec89[imo, ijr, :, :] = sa_89[ijr, :, :]
+        ratio89[imo, ijr, :, :] = ra_89[ijr, :, :]
+      
+    
+cumul = abs(grad_rat) + abs(spec23) + abs(spec89)
+cumul_month = np.zeros([M, ny, nx], float)
 for imo in range (0, M):
-    print 'month ' + month[imo]
-    ##############
-    # emissivity #
-    ##############
-    print 'read file for emiss'
-    fichier_emis = Dataset('/net/argos/data/parvati/lahlod/ARCTIC/monthly_GLACE/gridded_data/cartesian_grid/res_40/cartesian_grid_monthly_data_lamb_spec_near_nadir_AMSUA' + str(frequ) + '_' + month[imo] + '2009.nc', 'r', format='NETCDF3_CLASSIC')
-    xdist = fichier_emis.variables['longitude'][:]
-    ydist = fichier_emis.variables['latitude'][:]
-    day = fichier_emis.variables['days'][:]
-    emis_spec = fichier_emis.variables['e_spec'][:]
-    emis_lamb = fichier_emis.variables['e_lamb'][:]
-    fichier_emis.close()
-    print 'calculation of monthly mean'
-    for ilon in range (0, nx):
-        for ilat in range (0, ny):
-            emis_spec_moy[ilat, ilon, imo] = mean(emis_spec[ilat, ilon, :][nonzero(isnan(emis_spec[ilat, ilon, :]) == False)])
-            emis_lamb_moy[ilat, ilon, imo] = mean(emis_lamb[ilat, ilon, :][nonzero(isnan(emis_lamb[ilat, ilon, :]) == False)])
-    print 'calculation of monthly mean difference lamb-spec'
-    emis_diff[:, :, imo] = emis_lamb_moy[:, :, imo] - emis_spec_moy[:, :, imo]
-    #########
-    # ratio #
-    #########
-    print 'read file for ratio'
-    fichier_ratio = Dataset('/net/argos/data/parvati/lahlod/ARCTIC/monthly_GLACE/gridded_data/cartesian_grid/res_40/cartesian_grid_monthly_lamb-spec_ratio_near_nadir_AMSUA' + str(frequ) + '_' + month[imo] + '2009.nc', 'r', format='NETCDF3_CLASSIC')
-    xdist = fichier_ratio.variables['longitude'][:]
-    ydist = fichier_ratio.variables['latitude'][:]
-    emis_ratio[:, :, imo] = fichier_ratio.variables['emis_ratio'][:]
-    fichier_ratio.close()
+    for ilat in range (0, ny):
+        for ilon in range (0, nx):
+            cumul_month[imo, ilat, ilon] = mean(cumul[imo, :, ilat, ilon][nonzero(isnan(cumul[imo, :, ilat, ilon]) == False)])
+    
 
 
 
-###########################################
-# emissivity distribution after filtering #
-###########################################
+ion()
+x_ind, y_ind, z_ind, volume = arctic_map.arctic_map_lat50()
+x_coast = x_ind
+y_coast = y_ind
+z_coast = z_ind
+for imo in range (0, M):
+    map_cartesian_grid.draw_map_cartes_l(x_coast, y_coast, z_coast, volume, xvec, yvec, cumul_month[imo, :, :], 0., 0.9, 0.01, cm.s3pcpn_l_r, 'cumul')
+    title(month[imo] + ' 2009')
+    #plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA89_AMSUB89/maps/map_bias_AMSUA89-AMSUB89_arctic_' + month[imo] + '2009.png')
 
-hist_val_spec = np.zeros([50, M], float)
-hist_val_lamb = np.zeros([50, M], float)
-hist_val_ratio = np.zeros([50, M], float)
-hist_val_diff = np.zeros([50, M], float)
-corresp_val_spec = np.zeros([50, M], float)
-corresp_val_lamb = np.zeros([50, M], float)
-corresp_val_ratio = np.zeros([50, M], float)
-corresp_val_diff = np.zeros([50, M], float)
-emis_spec_f = np.zeros([7000, M], float)
-emis_lamb_f = np.zeros([7000, M], float)
-emis_ratio_f = np.zeros([7000, M], float)
-emis_diff_f = np.zeros([7000, M], float)
-X = np.zeros([7000, M], float)
-Y = np.zeros([7000, M], float)
-L_spec = np.zeros([M], int)
-for imo in range (0, M):
-    print 'month ' + month[imo]
-    # choice of spec emissivity as the threshold for the study // definition of x and y index corresponding to the points which emissivity value is over threshold
-    y_index_spec = np.where(emis_spec_moy[:, :, imo] >= lim_spec[imo])[0] 
-    Y[0 : len(y_index_spec), imo] = y_index_spec[:]
-    x_index_spec = np.where(emis_spec_moy[:, :, imo] >= lim_spec[imo])[1]
-    X[0 : len(x_index_spec), imo] = x_index_spec[:]
-    L_spec[imo] = len(x_index_spec) # = len(y_index)
-#print 'length of x and y index ', L_spec[imo]
-# definition of filtered values (vectors) with the previous threshold // values of filtered emissivity SPEC, LAMB, rate and difference LAMB-SPEC 
-for ii in range (0, L_spec[imo]):
-emis_spec_f[ii, imo] = emis_spec_moy[y_index_spec[ii], x_index_spec[ii], imo]
-emis_lamb_f[ii, imo] = emis_lamb_moy[y_index_spec[ii], x_index_spec[ii], imo]
-emis_ratio_f[ii, imo] = emis_ratio[y_index_spec[ii], x_index_spec[ii], imo]
-emis_diff_f[ii, imo] = emis_diff[y_index_spec[ii], x_index_spec[ii], imo]
-# definition of their distribution within the new filtered values
-hist_val_spec[:, imo] = hist(emis_spec_f, bins = 50, normed = True, histtype='step')[0]
-hist_val_lamb[:, imo] = hist(emis_lamb_f, bins = 50, normed = True, histtype='step')[0]
-hist_val_ratio[:, imo] = hist(emis_ratio_f, bins = 50, normed = True, histtype='step')[0]
-hist_val_diff[:, imo] = hist(emis_diff_f, bins = 50, normed = True, histtype='step')[0]
-for ibin in range (0, 50):
-corresp_val_spec[ibin, imo] = mean(hist(emis_spec_f, bins = 50, normed = True, histtype='step')[1][ibin : ibin + 2])
-corresp_val_lamb[ibin, imo] = mean(hist(emis_lamb_f, bins = 50, normed = True, histtype='step')[1][ibin : ibin + 2])
-corresp_val_ratio[ibin, imo] = mean(hist(emis_ratio_f, bins = 50, normed = True, histtype='step')[1][ibin : ibin + 2])
-corresp_val_diff[ibin, imo] = mean(hist(emis_diff_f, bins = 50, normed = True, histtype='step')[1][ibin : ibin + 2])
-
-
-return(emis_spec_f, emis_lamb_f, emis_ratio_f, emis_diff_f, L_spec)  
-
-
-
-'''
-########
-# plot #
-########
-c = np.array(['r', 'b', 'c', 'm', 'y', 'g'])
-figure()
-for imo in range (0, 6):
-plot(corresp_val[:, imo], hist_val[:, imo], c = str(c[imo]), label = str(month[imo]))
-
-grid()
-xlim(corresp_val.min() - 0.02, corresp_val.max() + 0.02)
-xlabel('emissivity parameter')
-ylabel('frequency of occurence')
-fontP = FontProperties()
-fontP.set_size('small')
-legend(prop = fontP)
-#plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/emiss_rate_AMSUA30_JANUARY-JUNE_2009.png')
-## plot six following months of spec emissivity histograms ##
-figure()
-for imo in range (6, M):
-plot(corresp_val[:, imo], hist_val[:, imo], c = str(c[imo - 6]), label = str(month[imo]))
-
-grid()
-xlim(corresp_val.min() - 0.02, corresp_val.max() + 0.02)
-xlabel('emissivity parameter')
-ylabel('frequency of occurence')
-legend(loc = 1, prop = fontP)
-#plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/emiss_rate_AMSUA30_JULY-DECEMBER_2009.png')
-'''