| 1 | #!/usr/bin/env python |
|---|
| 2 | # -*- coding: utf-8 -*- |
|---|
| 3 | import string |
|---|
| 4 | import numpy as np |
|---|
| 5 | import matplotlib.pyplot as plt |
|---|
| 6 | from pylab import * |
|---|
| 7 | from mpl_toolkits.basemap import Basemap |
|---|
| 8 | from mpl_toolkits.basemap import shiftgrid, cm |
|---|
| 9 | from netCDF4 import Dataset |
|---|
| 10 | import arctic_map # function to regrid coast limits |
|---|
| 11 | import cartesian_grid_test # function to convert grid from polar to cartesian |
|---|
| 12 | import scipy.special |
|---|
| 13 | import ffgrid2 |
|---|
| 14 | import map_ffgrid |
|---|
| 15 | from matplotlib import colors |
|---|
| 16 | from matplotlib.font_manager import FontProperties |
|---|
| 17 | import map_cartesian_grid |
|---|
| 18 | |
|---|
| 19 | |
|---|
| 20 | ############################### |
|---|
| 21 | # time period characteristics # |
|---|
| 22 | ############################### |
|---|
| 23 | MONTH = np.array(['01', '02', '03', '04', '05', '06', '07', '08', '09', '10', '11', '12']) |
|---|
| 24 | month = np.array(['JANUARY', 'FEBRUARY', 'MARCH', 'APRIL', 'MAY', 'JUNE', 'JULY', 'AUGUST', 'SEPTEMBER', 'OCTOBER', 'NOVEMBER', 'DECEMBER']) |
|---|
| 25 | month_day = np.array([31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31]) |
|---|
| 26 | M = len(month) |
|---|
| 27 | |
|---|
| 28 | |
|---|
| 29 | ######################## |
|---|
| 30 | # grid characteristics # |
|---|
| 31 | ######################## |
|---|
| 32 | x0 = -3000. # min limit of grid |
|---|
| 33 | x1 = 2500. # max limit of grid |
|---|
| 34 | dx = 40. |
|---|
| 35 | xvec = np.arange(x0, x1+dx, dx) |
|---|
| 36 | nx = len(xvec) |
|---|
| 37 | y0 = -3000. # min limit of grid |
|---|
| 38 | y1 = 3000. # max limit of grid |
|---|
| 39 | dy = 40. |
|---|
| 40 | yvec = np.arange(y0, y1+dy, dy) |
|---|
| 41 | ny = len(yvec) |
|---|
| 42 | |
|---|
| 43 | |
|---|
| 44 | |
|---|
| 45 | |
|---|
| 46 | ######################################################################## |
|---|
| 47 | # open NETCDF files of anomaly data and ice spec, lamb, diff and ratio # |
|---|
| 48 | ######################################################################## |
|---|
| 49 | frequ = np.array(['23', '30', '50', '89']) |
|---|
| 50 | |
|---|
| 51 | anom_spec = np.zeros([4, M, ny, nx], float) |
|---|
| 52 | anom_lamb = np.zeros([4, M, ny, nx], float) |
|---|
| 53 | anom_diff = np.zeros([4, M, ny, nx], float) |
|---|
| 54 | anom_ratio = np.zeros([4, M, ny, nx], float) |
|---|
| 55 | |
|---|
| 56 | ice_spec = np.zeros([4, M, ny, nx], float) |
|---|
| 57 | ice_lamb = np.zeros([4, M, ny, nx], float) |
|---|
| 58 | ice_diff = np.zeros([4, M, ny, nx], float) |
|---|
| 59 | ice_ratio = np.zeros([4, M, ny, nx], float) |
|---|
| 60 | |
|---|
| 61 | for ifr in range (0, 4): |
|---|
| 62 | for imo in range (0, M): |
|---|
| 63 | print 'stack in file month ' + str(month[imo]) |
|---|
| 64 | print 'open anomaly file' |
|---|
| 65 | fichier_anom = Dataset('/net/argos/data/parvati/lahlod/ARCTIC/AMSUA_ice_class/sub_classification/cartesian_grid_map_sea_ice_extent_with-AMSUA23-and-30_param_anomaly_' + month[imo] + '2009_AMSUA' + str(frequ[ifr]) + '.nc', 'r', format='NETCDF3_CLASSIC') |
|---|
| 66 | anom_spec[ifr, imo, :, :] = fichier_anom.variables['spec_anomaly'][:] |
|---|
| 67 | anom_lamb[ifr, imo, :, :] = fichier_anom.variables['lamb_anomay'][:] |
|---|
| 68 | anom_diff[ifr, imo, :, :] = fichier_anom.variables['diff_anomaly'][:] |
|---|
| 69 | anom_ratio[ifr, imo, :, :] = fichier_anom.variables['ratio_anomaly'][:] |
|---|
| 70 | fichier_anom.close() |
|---|
| 71 | print 'open ice file' |
|---|
| 72 | fichier_ice = Dataset('/net/argos/data/parvati/lahlod/ARCTIC/AMSUA_ice_class/sub_classification/cartesian_grid_map_sea_ice_extent_with-AMSUA23-and-30_' + month[imo] + '2009_AMSUA' + str(frequ[ifr]) + '_spec_lamb_thresholds.nc', 'r', format='NETCDF3_CLASSIC') |
|---|
| 73 | ice_spec[ifr, imo, :, :] = fichier_ice.variables['spec_ice_area'][:] |
|---|
| 74 | ice_lamb[ifr, imo, :, :] = fichier_ice.variables['lamb_ice_area'][:] |
|---|
| 75 | ice_diff[ifr, imo, :, :] = fichier_ice.variables['diff_ice_area'][:] |
|---|
| 76 | ice_ratio[ifr, imo, :, :] = fichier_ice.variables['ratio_ice_area'][:] |
|---|
| 77 | fichier_ice.close() |
|---|
| 78 | |
|---|
| 79 | |
|---|
| 80 | ''' |
|---|
| 81 | ############### |
|---|
| 82 | # map anomaly # |
|---|
| 83 | ############### |
|---|
| 84 | print 'map monthly anomaly' |
|---|
| 85 | ion() |
|---|
| 86 | x_ind, y_ind, z_ind, volume = arctic_map.arctic_map_lat50() |
|---|
| 87 | x_coast = x_ind |
|---|
| 88 | y_coast = y_ind |
|---|
| 89 | z_coast = z_ind |
|---|
| 90 | ifr = 0 # anomaly of emissivity ratio from AMSUA89GHz database |
|---|
| 91 | for imo in range (0, M): |
|---|
| 92 | print 'month ' + str(month[imo]) |
|---|
| 93 | map_cartesian_grid.draw_map_cartes_l(x_coast, y_coast, z_coast, volume, xvec, yvec, anom_spec[ifr, imo, :, :], anom_spec[ifr, imo,:,:][nonzero(isnan(anom_spec[ifr, imo,:,:]) == False)].min(), anom_spec[ifr, imo, :, :][nonzero(isnan(anom_spec[ifr, imo , :, :]) == False)].max(), 0.001, cm.s3pcpn_l_r, 'Anomaly of emissivity spec') |
|---|
| 94 | title('AMSUA ' + str(frequ[ifr]) + 'GHz - ' + str(month[imo]) + ' 2009') |
|---|
| 95 | print 'save figure in .png file' |
|---|
| 96 | plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/sub_classification/maps_sea_ice_extent/emiss_anomaly/spec/anomaly_spec_map_AMSUA' + str(frequ[ifr]) + '_with_AMSUA23-and-30-spec_'+str(month[imo])+'_2009.png') |
|---|
| 97 | #cm.s3pcpn_l_r |
|---|
| 98 | ''' |
|---|
| 99 | |
|---|
| 100 | |
|---|
| 101 | ################################################ |
|---|
| 102 | # calculate gradient ratio for two frequencies # |
|---|
| 103 | ################################################ |
|---|
| 104 | print 'calculate gradient ratio' |
|---|
| 105 | grad_ratio = np.zeros([M, ny, nx], float) |
|---|
| 106 | hist_ratio = np.zeros([M, 50], float) |
|---|
| 107 | corresp_ratio = np.zeros([M, 50], float) |
|---|
| 108 | for imo in range (0, M): |
|---|
| 109 | print 'month ' + month[imo] |
|---|
| 110 | for ilat in range (0, ny): |
|---|
| 111 | for ilon in range (0, nx): |
|---|
| 112 | grad_ratio[imo, ilat, ilon] = (ice_spec[0, imo, ilat, ilon] - ice_spec[3, imo, ilat, ilon]) / (ice_spec[0, imo, ilat, ilon] + ice_spec[3, imo, ilat, ilon]) |
|---|
| 113 | ''' |
|---|
| 114 | # compute histogram distribution of gradient ratio to find characteristic threshold |
|---|
| 115 | grad_ratio_vec = reshape(grad_ratio[imo, :, :], size(grad_ratio[imo, :, :]))[nonzero(isnan(reshape(grad_ratio[imo, :, :], size(grad_ratio[imo, :, :]))) == False)] |
|---|
| 116 | hist_ratio[imo, :] = hist(grad_ratio_vec, bins = 50, normed = True, histtype='step')[0] |
|---|
| 117 | for ibin in range (0, 50): |
|---|
| 118 | corresp_ratio[imo, ibin] = mean(hist(grad_ratio_vec, bins = 50, normed = True, histtype='step')[1][ibin : ibin + 2]) |
|---|
| 119 | |
|---|
| 120 | # plot histogram |
|---|
| 121 | c = np.array(['r', 'b', 'c', 'm', 'y', 'g']) |
|---|
| 122 | figure() |
|---|
| 123 | for imo in range (0, 6): |
|---|
| 124 | plot(corresp_ratio[imo], hist_ratio[imo], c = str(c[imo]), label = str(month[imo])) |
|---|
| 125 | |
|---|
| 126 | grid() |
|---|
| 127 | xlim(corresp_ratio.min(), corresp_ratio.max()) |
|---|
| 128 | xlabel('emissivity ratio') |
|---|
| 129 | ylabel('frequency of occurence') |
|---|
| 130 | fontP = FontProperties() |
|---|
| 131 | fontP.set_size('small') |
|---|
| 132 | legend(loc = 2, prop = fontP) |
|---|
| 133 | #plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/emiss_spec_AMSUA' + str(frequ) + '_JANUARY-JUNE_2009.png') |
|---|
| 134 | ## plot six following months of spec emissivity histograms ## |
|---|
| 135 | figure() |
|---|
| 136 | for imo in range (6, M): |
|---|
| 137 | plot(corresp_ratio[imo], hist_ratio[imo], c = str(c[imo - 6]), label = str(month[imo])) |
|---|
| 138 | |
|---|
| 139 | grid() |
|---|
| 140 | xlim(corresp_ratio.min(), corresp_ratio.max()) |
|---|
| 141 | xlabel('emissivity ratio') |
|---|
| 142 | ylabel('frequency of occurence') |
|---|
| 143 | legend(loc = 1, prop = fontP) |
|---|
| 144 | #plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/emiss_spec_AMSUA' + str(frequ) + '_JULY-DECEMBER_2009.png') |
|---|
| 145 | ''' |
|---|
| 146 | |
|---|
| 147 | ################################ |
|---|
| 148 | # map cumulation of parameters # |
|---|
| 149 | ################################ |
|---|
| 150 | cumul1 = anom_spec[0, :, :, :] - grad_ratio # correlation |
|---|
| 151 | cumul2 = anom_spec[3, :, :, :] + grad_ratio # anti correlation |
|---|
| 152 | cumul3 = anom_spec[3, :, :, :] + anom_ratio[3, :, :, :] # anti correlation |
|---|
| 153 | cumul4 = abs(anom_spec[0, :, :, :]) + abs(anom_spec[3, :, :, :]) + abs(grad_ratio) + abs(anom_ratio[3, :, :, :]) |
|---|
| 154 | |
|---|
| 155 | |
|---|
| 156 | |
|---|
| 157 | |
|---|
| 158 | ion() |
|---|
| 159 | x_ind, y_ind, z_ind, volume = arctic_map.arctic_map_lat50() |
|---|
| 160 | x_coast = x_ind |
|---|
| 161 | y_coast = y_ind |
|---|
| 162 | z_coast = z_ind |
|---|
| 163 | for imo in range (0, M): |
|---|
| 164 | ''' |
|---|
| 165 | # map cumul1 |
|---|
| 166 | map_cartesian_grid.draw_map_cartes_l(x_coast, y_coast, z_coast, volume, xvec, yvec, cumul1[imo, :, :], -0.25, 0.14, 0.01, cm.s3pcpn_l_r, 'Cumulation [emis spec anomaly 23 - grad ratio 23-89]') |
|---|
| 167 | title('AMSUA - ' + month[imo] + ' 2009') |
|---|
| 168 | plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/sub_classification/maps_sea_ice_extent/cumul_params/map_cumulation_spec-anom-23_grad-ratio_' + month[imo] +'2009.png') |
|---|
| 169 | # map cumul2 |
|---|
| 170 | map_cartesian_grid.draw_map_cartes_l(x_coast, y_coast, z_coast, volume, xvec, yvec, cumul2[imo, :, :], -0.17, 0.17, 0.01, cm.s3pcpn_l_r, 'Cumulation [emis spec anomaly 89 + grad ratio 23-89]') |
|---|
| 171 | title('AMSUA - ' + month[imo] + ' 2009') |
|---|
| 172 | plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/sub_classification/maps_sea_ice_extent/cumul_params/map_cumulation_spec-anom-89_grad-ratio_' + month[imo] +'2009.png') |
|---|
| 173 | # map cumul3 |
|---|
| 174 | map_cartesian_grid.draw_map_cartes_l(x_coast, y_coast, z_coast, volume, xvec, yvec, cumul3[imo, :, :], -5.46, 2.14, 0.01, cm.s3pcpn_l_r, 'Cumulation [emis spec anomaly 89 + emis ratio 89]') |
|---|
| 175 | title('AMSUA - ' + month[imo] + ' 2009') |
|---|
| 176 | plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/sub_classification/maps_sea_ice_extent/cumul_params/map_cumulation_spec-anom-89_ratio-anom-89_' + month[imo] +'2009.png') |
|---|
| 177 | ''' |
|---|
| 178 | # map cumul4 |
|---|
| 179 | map_cartesian_grid.draw_map_cartes_l(x_coast, y_coast, z_coast, volume, xvec, yvec, cumul4[imo, :, :], cumul4[imo, :, :][nonzero(isnan(cumul4[imo, :, :]) == False)].min(), cumul4[imo, :, :][nonzero(isnan(cumul4[imo, :, :]) == False)].max(), 0.1, cm.s3pcpn_l_r, 'Cumulation all parameters') |
|---|
| 180 | title('AMSUA - ' + month[imo] + ' 2009') |
|---|
| 181 | plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/sub_classification/maps_sea_ice_extent/cumul_params/map_cumulation_all_parameters_' + month[imo] +'2009.png') |
|---|
| 182 | |
|---|
| 183 | |
|---|
| 184 | |
|---|
| 185 | |
|---|
| 186 | |
|---|
| 187 | ''' |
|---|
| 188 | figure() |
|---|
| 189 | scatter(x_coast, y_coast, c = z_coast, s = volume) |
|---|
| 190 | #cmap = cm.s3pcpn_l_r |
|---|
| 191 | levels = np.arange(-1, 11, 1) |
|---|
| 192 | cs = contour(xvec, yvec, cumul[imo, :, :]) |
|---|
| 193 | cs = contour(xvec, yvec, new_anom_spec0[imo, :, :], levels) |
|---|
| 194 | cs = contour(xvec, yvec, new_anom_ratio[imo, :, :], levels) |
|---|
| 195 | cs = contour(xvec, yvec, new_grad_ratio[imo, :, :], levels) |
|---|
| 196 | cbar = colorbar(cs, ticks = [-1, 0, 1, 2, 3, 4, 5]) |
|---|
| 197 | cbar.ax.set_yticklabels(['emis_spec_89', 'emis_spec_23', 'emis_ratio_89', 'grad_ratio_23-89', '-']) |
|---|
| 198 | #cbar.set_label(text) |
|---|
| 199 | xlim(-3500, 2700.) |
|---|
| 200 | ylim(-4000, 2800.) |
|---|
| 201 | ''' |
|---|
| 202 | |
|---|
| 203 | |
|---|
| 204 | |
|---|
| 205 | |
|---|
| 206 | |
|---|
| 207 | |
|---|
| 208 | |
|---|
| 209 | ''' |
|---|
| 210 | ############################################ |
|---|
| 211 | # time evolution (monthly) in a given zone # |
|---|
| 212 | ############################################ |
|---|
| 213 | # read monthly std of emis |
|---|
| 214 | ifr = 3 |
|---|
| 215 | 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') |
|---|
| 216 | |
|---|
| 217 | |
|---|
| 218 | |
|---|
| 219 | #zone 1 : yi = 400. // yf = 720. // xi = -1080. // xf = -920. (Canadian Archipelago) |
|---|
| 220 | #zone 2 : yi = 920. // yf = 1280. // xi = -720. // xf = -560. (Chukchi Sea, Bering Detroit) |
|---|
| 221 | # select borders of zone |
|---|
| 222 | yi = 920. |
|---|
| 223 | yf = 1280. |
|---|
| 224 | xi = -720. |
|---|
| 225 | xf = -560. |
|---|
| 226 | #find corresponding index in xvec and yvec |
|---|
| 227 | xxi = np.where(xvec == xi)[0][0] |
|---|
| 228 | xxf = np.where(xvec == xf)[0][0] |
|---|
| 229 | yyi = np.where(yvec == yi)[0][0] |
|---|
| 230 | yyf = np.where(yvec == yf)[0][0] |
|---|
| 231 | |
|---|
| 232 | # define vectors |
|---|
| 233 | anom_spec0_zone = np.zeros([M], float) |
|---|
| 234 | anom_spec3_zone = np.zeros([M], float) |
|---|
| 235 | anom_ratio_zone = np.zeros([M], float) |
|---|
| 236 | grad_ratio_zone = np.zeros([M], float) |
|---|
| 237 | std_as0 = np.zeros([M], float) |
|---|
| 238 | std_as3 = np.zeros([M], float) |
|---|
| 239 | std_ar = np.zeros([M], float) |
|---|
| 240 | std_gr = np.zeros([M], float) |
|---|
| 241 | for imo in range (0, M): |
|---|
| 242 | # anom spec emis from 23GHz |
|---|
| 243 | anom_spec0_vec = reshape(anom_spec[0, imo, yyi : yyf + 1, xxi : xxf + 1], size(anom_spec[0, imo, yyi : yyf + 1, xxi : xxf + 1])) |
|---|
| 244 | if (len(anom_spec0_vec[nonzero(isnan(anom_spec0_vec) == True)]) != 0): |
|---|
| 245 | anom_spec0_zone[imo] = NaN |
|---|
| 246 | std_as0[imo] = NaN |
|---|
| 247 | anom_spec3_zone[imo] = NaN |
|---|
| 248 | std_as3[imo] = NaN |
|---|
| 249 | anom_ratio_zone[imo] = NaN |
|---|
| 250 | std_ar[imo] = NaN |
|---|
| 251 | grad_ratio_zone[imo] = NaN |
|---|
| 252 | std_gr[imo] = NaN |
|---|
| 253 | continue |
|---|
| 254 | else: |
|---|
| 255 | anom_spec0_zone[imo] = anom_spec0_vec.mean() |
|---|
| 256 | std_as0[imo] = sqrt((1. / (size(anom_spec[0, imo, yyi : yyf + 1, xxi : xxf + 1]) - 1.)) * sum((anom_spec0_vec[:] - anom_spec0_zone[imo])**2)) |
|---|
| 257 | # anom spec emis from 89GHz |
|---|
| 258 | anom_spec3_vec = reshape(anom_spec[3, imo, yyi : yyf + 1, xxi : xxf + 1], size(anom_spec[3, imo, yyi : yyf + 1, xxi : xxf + 1])) |
|---|
| 259 | anom_spec3_zone[imo] = anom_spec3_vec.mean() |
|---|
| 260 | std_as3[imo] = sqrt((1. / (size(anom_spec[3, imo, yyi : yyf + 1, xxi : xxf + 1]) - 1.)) * sum((anom_spec3_vec[:] - anom_spec3_zone[imo])**2)) |
|---|
| 261 | # anom emis ratio from 89GHz |
|---|
| 262 | anom_ratio_vec = reshape(anom_ratio[3, imo, yyi : yyf + 1, xxi : xxf + 1], size(anom_ratio[3, imo, yyi : yyf + 1, xxi : xxf + 1])) |
|---|
| 263 | anom_ratio_zone[imo] = anom_ratio_vec.mean() |
|---|
| 264 | std_ar[imo] = sqrt((1. / (size(anom_ratio[3, imo, yyi : yyf + 1, xxi : xxf + 1]) - 1.)) * sum((anom_ratio_vec[:] - anom_ratio_zone[imo])**2)) |
|---|
| 265 | # gradient ratio from 23 and 89GHz |
|---|
| 266 | grad_ratio_vec = reshape(grad_ratio[imo, yyi : yyf + 1, xxi : xxf + 1], size(grad_ratio[imo, yyi : yyf + 1, xxi : xxf + 1])) |
|---|
| 267 | grad_ratio_zone[imo] = grad_ratio_vec.mean() |
|---|
| 268 | std_gr[imo] = sqrt((1. / (size(grad_ratio[imo, yyi : yyf + 1, xxi : xxf + 1]) - 1.)) * sum((grad_ratio_vec[:] - grad_ratio_zone[imo])**2)) |
|---|
| 269 | |
|---|
| 270 | |
|---|
| 271 | figure() |
|---|
| 272 | fontP = FontProperties() |
|---|
| 273 | fontP.set_size('small') |
|---|
| 274 | subplot(2,1,1) |
|---|
| 275 | errorbar(arange(0, M, 1), anom_spec0_zone, yerr = std_as0, label = 'emis spec anomaly 23GHz') |
|---|
| 276 | errorbar(arange(0, M, 1), anom_spec3_zone, yerr = std_as3, label = 'emis spec anomaly 89GHz') |
|---|
| 277 | errorbar(arange(0, M, 1), grad_ratio_zone, yerr = std_gr, label = 'gradient ratio 23-89GHz') |
|---|
| 278 | plot(arange(0, M, 1), np.zeros([M], float), '--k') |
|---|
| 279 | legend(loc = 3, prop = fontP) |
|---|
| 280 | grid() |
|---|
| 281 | xticks(arange(0, M, 1), month, rotation = 15) |
|---|
| 282 | xlim(-1, M) |
|---|
| 283 | subplot(2,1,2) |
|---|
| 284 | errorbar(arange(0, M, 1), anom_ratio_zone, yerr = std_ar, label = 'emis ratio anomaly 89GHz') |
|---|
| 285 | plot(arange(0, M, 1), np.zeros([M], float), '--k') |
|---|
| 286 | legend(loc = 2, prop = fontP) |
|---|
| 287 | grid() |
|---|
| 288 | xticks(arange(0, M, 1), month, rotation = 15) |
|---|
| 289 | xlim(-1, M) |
|---|
| 290 | plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/sub_classification/study_by_zones/monthly_evolution_emis_params_zone_North_Beaufort_Sea.png') |
|---|
| 291 | ### map the selected zone ### |
|---|
| 292 | map_cartesian_grid.draw_map_cartes_l(x_coast, y_coast, z_coast, volume, xvec[xxi : xxf], yvec[yyi : yyf], grad_ratio[imo, yyi : yyf, xxi : xxf], grad_ratio[imo, yyi : yyf, xxi : xxf][nonzero(isnan(grad_ratio[imo, yyi : yyf, xxi : xxf]) == False)].min(), grad_ratio[imo, yyi : yyf, xxi : xxf][nonzero(isnan(grad_ratio[imo, yyi : yyf, xxi : xxf]) == False)].max(), 0.001, cm.s3pcpn_l_r, 'selected zone') |
|---|
| 293 | plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/sub_classification/study_by_zones/map_zone_North_Beaufort_Sea.png') |
|---|
| 294 | ''' |
|---|
| 295 | |
|---|
| 296 | |
|---|
| 297 | |
|---|
| 298 | ''' |
|---|
| 299 | ############################################################################ |
|---|
| 300 | # calculate correlation between anom ratio and gradient ratio (spec emiss) # |
|---|
| 301 | ############################################################################ |
|---|
| 302 | corr_mat = np.zeros([M, 4, 4], float) |
|---|
| 303 | for imo in range (0, M): |
|---|
| 304 | a = reshape(anom_spec[0, imo, :, :], size(anom_spec[0, imo, :, :]))[nonzero(isnan(reshape(anom_spec[0, imo, :, :], size(anom_spec[0, imo, :, :]))) == False)] |
|---|
| 305 | b = reshape(anom_spec[3, imo, :, :], size(anom_spec[3, imo, :, :]))[nonzero(isnan(reshape(anom_spec[3, imo, :, :], size(anom_spec[3, imo, :, :]))) == False)] |
|---|
| 306 | c = reshape(anom_ratio[3, imo, :, :], size(anom_ratio[3, imo, :, :]))[nonzero(isnan(reshape(anom_ratio[3, imo, :, :], size(anom_ratio[3, imo, :, :]))) == False)] |
|---|
| 307 | d = reshape(grad_ratio[imo, :, :], size(grad_ratio[imo, :, :]))[nonzero(isnan(reshape(grad_ratio[imo, :, :], size(grad_ratio[imo, :, :]))) == False)] |
|---|
| 308 | params = np.array([a, b, c, d]) |
|---|
| 309 | for ii in range (0, 4): |
|---|
| 310 | for jj in range (0, 4): |
|---|
| 311 | corr_mat[imo, ii, jj] = corrcoef(params[ii], params[jj])[0][1] |
|---|
| 312 | |
|---|
| 313 | figure() |
|---|
| 314 | fontP = FontProperties() |
|---|
| 315 | fontP.set_size('small') |
|---|
| 316 | plot(arange(0, M, 1), corr_mat[:, 0, 1], 'r', label = 'spec anomaly 23GHz - spec anomaly 89GHz') |
|---|
| 317 | plot(arange(0, M, 1), corr_mat[:, 0, 2], 'g', label = 'spec anomaly 23GHz - ratio anomaly 89GHz') |
|---|
| 318 | plot(arange(0, M, 1), corr_mat[:, 0, 3], 'b', label = 'spec anomaly 23GHz - gradient ratio') |
|---|
| 319 | plot(arange(0, M, 1), corr_mat[:, 1, 2], 'c', label = 'spec anomaly 89GHz - ratio anomaly 89GHz') |
|---|
| 320 | plot(arange(0, M, 1), corr_mat[:, 1, 3], 'm', label = 'spec anomaly 89GHz - gradient ratio') |
|---|
| 321 | plot(arange(0, M, 1), corr_mat[:, 2, 3], 'y', label = 'ratio anomaly 89GHz - gradient ratio') |
|---|
| 322 | plot(arange(0, M, 1), zeros([M], float), '--k') |
|---|
| 323 | xlim(-1, M) |
|---|
| 324 | ylim(-2, 1.) |
|---|
| 325 | legend(loc = 8, prop = fontP) |
|---|
| 326 | xticks(arange(0, M, 1), month, rotation = 25) |
|---|
| 327 | yticks(np.arange(-2., 1., 0.1)) |
|---|
| 328 | grid() |
|---|
| 329 | plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/sub_classification/correlation_params_2009.png') |
|---|
| 330 | ''' |
|---|
| 331 | |
|---|
| 332 | ''' |
|---|
| 333 | corr_mat = np.zeros([4, ny, nx], float) |
|---|
| 334 | for ifr in range (0, 4): |
|---|
| 335 | for ilat in range (0, ny): |
|---|
| 336 | for ilon in range (0, nx): |
|---|
| 337 | corr_mat[ifr, ilat, ilon] = corrcoef(grad_ratio[:, ilat, ilon], anom_spec[ifr, :, ilat, ilon])[0][1] |
|---|
| 338 | |
|---|
| 339 | for ifr in range (0, 4): |
|---|
| 340 | map_cartesian_grid.draw_map_cartes_l(x_coast, y_coast, z_coast, volume, xvec, yvec, corr_mat[ifr, :, :], -1., 1., 0.01, cm.s3pcpn_l_r, 'correlation [emis ratio anomaly - gradient ratio]') |
|---|
| 341 | title('AMSUA ' + str(frequ[ifr]) + 'GHz - year 2009') |
|---|
| 342 | plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/sub_classification/maps_sea_ice_extent/correlation_map_grad-ratio_emis-spec-anom_AMSUA' + str(frequ[ifr]) + '_with_AMSUA23-and-30-spec_year2009.png') |
|---|
| 343 | ''' |
|---|
| 344 | |
|---|
| 345 | |
|---|
| 346 | |
|---|
| 347 | |
|---|
| 348 | |
|---|
| 349 | |
|---|
| 350 | |
|---|
| 351 | |
|---|
![(please configure the [header_logo] section in trac.ini)](/twice/chrome/site/your_project_logo.png){kind=logo}
