![(please configure the [header_logo] section in trac.ini)](/twice/chrome/site/your_project_logo.png){kind=logo}
- Timestamp:
- 07/30/14 18:21:31 (10 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/src/scripts_Laura/ARCTIC/Travail_CEN/ice_delimit_AMSUB_data.py
r46 r47 41 41 42 42 43 frequ = 89 43 44 44 45 ##################### … … 56 57 ### emissivity ### 57 58 print 'open emissivity file' 58 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_AMSU A30_' + str(month[imo]) + '2009.nc', 'r', format='NETCDF3_CLASSIC')59 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_AMSUB' + str(frequ) + '_' + str(month[imo]) + '2009.nc', 'r', format='NETCDF3_CLASSIC') 59 60 xdist = fichier_emis.variables['longitude'][:] 60 61 ydist = fichier_emis.variables['latitude'][:] … … 70 71 ### monthly emissivity ratio ### 71 72 print 'open emissivity ratio file' 72 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_AMSU A30_' + str(month[imo]) + '2009.nc', 'r', format='NETCDF3_CLASSIC')73 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_AMSUB' + str(frequ) + '_' + str(month[imo]) + '2009.nc', 'r', format='NETCDF3_CLASSIC') 73 74 xdist = fichier_ratio.variables['longitude'][:] 74 75 ydist = fichier_ratio.variables['latitude'][:] … … 94 95 95 96 c = np.array(['r', 'b', 'c', 'm', 'y', 'g']) 96 limit_coef_spec = np.array([0.6, 0.6, 0.7, 0.8, 0.75]) # i1 = AMSUA23GHz / i2 = AMSUA30GHz / i3 = AMSUA50GHz / i4 = AMSUA89GHz / i5 = AMSUB89GHz 97 limit_coef_lamb = np.array([0.6, 0.6, 0.8, 0.8, 0.77]) # i1 = AMSUA23GHz / i2 = AMSUA30GHz / i3 = AMSUA50GHz / i4 = AMSUA89GHz / i5 = AMSUB89GHz 98 idata = 1 97 #limit_coef_spec = np.array([0.6, 0.6, 0.7, 0.8, 0.75]) # i1 = AMSUA23GHz / i2 = AMSUA30GHz / i3 = AMSUA50GHz / i4 = AMSUA89GHz / i5 = AMSUB89GHz 98 #limit_coef_lamb = np.array([0.6, 0.6, 0.8, 0.8, 0.77]) # i1 = AMSUA23GHz / i2 = AMSUA30GHz / i3 = AMSUA50GHz / i4 = AMSUA89GHz / i5 = AMSUB89GHz 99 #idata = 0 100 fontP = FontProperties() 101 fontP.set_size('small') 99 102 print 'distribution of emissivity values' 100 103 ################################### … … 119 122 xlabel('emissivity SPEC') 120 123 ylabel('frequency of occurence') 121 fontP = FontProperties()122 fontP.set_size('small')123 124 legend(loc = 2, prop = fontP) 124 plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/emiss_spec_AMSUA30_JANUARY-JUNE_2009.png')125 #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') 125 126 ## plot six following months of spec emissivity histograms ## 126 127 figure() … … 133 134 ylabel('frequency of occurence') 134 135 legend(loc = 9, prop = fontP) 135 plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/emiss_spec_AMSUA30_JULY-DECEMBER_2009.png')136 #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') 136 137 137 138 # find the emissivity value corresponding to the threshold of ice/no_ice limit … … 139 140 for imo in range (0, M): 140 141 print 'month ' + str(month[imo]) 141 bb = np.where((corresp_emis_spec[imo, :] < limit_coef_spec[idata]))[0] # only consider low emissivity values (open sea) // CHANGE EMISSIVITY VALUE ACCORDING TO FREQUENCY142 aa = np.where(hist_vals_spec[imo, bb] < max(hist_vals_spec[imo, bb]) / 2.)[0] # of these latter values, only consider emissivity values lower than 1/4*peak emissivity143 cc = np.where(hist_vals_spec[imo, bb] == max(hist_vals_spec[imo, bb]))[0] # which emissivity index corresponds to peak emissivity144 dd = np.where(aa > cc[0])[0]145 emis_lim_spec[imo] = corresp_emis_spec[imo, bb][aa][ dd[0]]142 bb = np.where((corresp_emis_spec[imo, :] > 0.7) & (corresp_emis_spec[imo, :] < 0.77))[0] # only consider low emissivity values (open sea) // CHANGE EMISSIVITY VALUE ACCORDING TO FREQUENCY 143 aa = np.where(hist_vals_spec[imo, bb] == min(hist_vals_spec[imo, bb]))[0] # of these latter values, only consider emissivity values lower than 1/4*peak emissivity 144 #cc = np.where(hist_vals_spec[imo, bb] == max(hist_vals_spec[imo, bb]))[0] # which emissivity index corresponds to peak emissivity 145 #dd = np.where(aa > cc[0])[0] 146 emis_lim_spec[imo] = corresp_emis_spec[imo, bb][aa][0] 146 147 147 148 # plot monthly evolution of emissivity spec threshold used for delimitation … … 153 154 legend(loc = 2) 154 155 xlim(-1, M) 155 plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSU A/emis_lim_frequ_SPEC_AMSUA30_2009.png')156 plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUB/emis_lim_frequ_SPEC_AMSUB' + str(frequ) + '_2009.png') 156 157 157 158 … … 177 178 ylabel('frequency of occurence') 178 179 legend(loc = 2, prop = fontP) 179 plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/emiss_lamb_AMSUA30_JANUARY-JUNE_2009.png')180 #plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/emiss_lamb_AMSUA' + str(frequ) + '_JANUARY-JUNE_2009.png') 180 181 ## plot six following months of spec emissivity histograms ## 181 182 figure() … … 188 189 ylabel('frequency of occurence') 189 190 legend(loc = 9, prop = fontP) 190 plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/emiss_lamb_AMSUA30_JULY-DECEMBER_2009.png')191 #plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/emiss_lamb_AMSUA' + str(frequ) + '_JULY-DECEMBER_2009.png') 191 192 192 193 # find the emissivity value corresponding to the threshold of ice/no_ice limit … … 194 195 for imo in range (0, M): 195 196 print 'month ' + str(month[imo]) 196 bb = np.where((corresp_emis_lamb[imo, :] < limit_coef_lamb[idata]))[0] # only consider low emissivity values (open sea) // CHANGE EMISSIVITY VALUE ACCORDING TO FREQUENCY197 aa = np.where(hist_vals_lamb[imo, bb] < max(hist_vals_lamb[imo, bb]) / 2.)[0] # of these latter values, only consider emissivity values lower than 1/4*peak emissivity198 cc = np.where(hist_vals_lamb[imo, bb] == max(hist_vals_lamb[imo, bb]))[0] # which emissivity index corresponds to peak emissivity199 dd = np.where(aa > cc[0])[0]200 emis_lim_lamb[imo] = corresp_emis_lamb[imo, bb][aa][ dd[0]]197 bb = np.where((corresp_emis_lamb[imo, :] > 0.72) & (corresp_emis_lamb[imo, :] < 0.8))[0] # only consider low emissivity values (open sea) // CHANGE EMISSIVITY VALUE ACCORDING TO FREQUENCY 198 aa = np.where(hist_vals_lamb[imo, bb] == min(hist_vals_lamb[imo, bb]))[0] # of these latter values, only consider emissivity values lower than 1/4*peak emissivity 199 #cc = np.where(hist_vals_lamb[imo, bb] == max(hist_vals_lamb[imo, bb]))[0] # which emissivity index corresponds to peak emissivity 200 #dd = np.where(aa > cc[0])[0] 201 emis_lim_lamb[imo] = corresp_emis_lamb[imo, bb][aa][0] 201 202 202 203 # plot monthly evolution of emissivity spec threshold used for delimitation … … 208 209 legend(loc = 2) 209 210 xlim(-1, M) 210 plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSU A/emis_lim_frequ_LAMB_AMSUA30_2009.png')211 plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUB/emis_lim_frequ_LAMB_AMSUB' + str(frequ) + '_2009.png') 211 212 212 213 … … 233 234 ylabel('frequency of occurence') 234 235 legend(prop = fontP) 235 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')236 #plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/emiss_rate_AMSUA'+str(frequ)+'_JANUARY-JUNE_2009.png') 236 237 ## plot six following months of spec emissivity histograms ## 237 238 figure() … … 244 245 ylabel('frequency of occurence') 245 246 legend(loc = 9, prop = fontP) 246 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')247 #plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/emiss_rate_AMSUA'+str(frequ)+'_JULY-DECEMBER_2009.png') 247 248 248 249 # no definition of threshold for ICE / NO_ICE delimitation because no apparent distinction signal // use of this signal for internal ice classification … … 280 281 legend(loc = 2, prop = fontP) 281 282 xlim(-1, M) 282 plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSU A/emis_lim_frequ_spec_and_lamb_AMSUA30_2009.png')283 plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUB/emis_lim_frequ_spec_and_lamb_AMSUB' + str(frequ) + '_2009.png') 283 284 284 285 … … 304 305 else: 305 306 if (emis_spec_month[imo, ilat, ilon] <= emis_lim_spec[imo]): # use the monthly SPEC emissivity threshold 306 ice_zone_spec[imo, ilat, ilon] = 0.307 ice_zone_spec[imo, ilat, ilon] = NaN 307 308 else: 308 ice_zone_spec[imo, ilat, ilon] = 1.309 map_cartesian_grid.draw_map_cartes_l(x_coast, y_coast, z_coast, volume, xdist, ydist, ice_zone_spec[ imo, :, :], -1, 2., 1., colors.ListedColormap(['0.5', 'blue']), 'Ice limit(threshold : emis_SPEC > ' + str(emis_lim_spec[imo])[0:6] + ')')310 title(str(month[imo]) + ' 2009 - AMSUA 30GHz')311 plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/AMSUA 30_ice_emis_spec_thresh' + '_' + month[imo] + '2009.png')309 ice_zone_spec[imo, ilat, ilon] = emis_spec_month[imo, ilat, ilon] 310 map_cartesian_grid.draw_map_cartes_l(x_coast, y_coast, z_coast, volume, xdist, ydist, ice_zone_spec[8, :, :], 0.5, 1.02, 0.01, cm.s3pcpn_l_r, 'Sea ice emissivity spec (threshold : emis_SPEC > ' + str(emis_lim_spec[imo])[0:6] + ')') 311 title(str(month[imo]) + ' 2009 - AMSUA ' + str(frequ) + 'GHz') 312 plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/AMSUA' + str(frequ) + '_ice_emis_spec_thresh' + '_' + month[imo] + '2009.png') 312 313 313 314 #### LAMB emiss #### … … 322 323 else: 323 324 if (emis_lamb_month[imo, ilat, ilon] <= emis_lim_lamb[imo]): # use the monthly SPEC emissivity threshold 324 ice_zone_lamb[imo, ilat, ilon] = 0.325 ice_zone_lamb[imo, ilat, ilon] = NaN 325 326 else: 326 ice_zone_lamb[imo, ilat, ilon] = 1.327 map_cartesian_grid.draw_map_cartes_l(x_coast, y_coast, z_coast, volume, xdist, ydist, ice_zone_lamb[ imo, :, :], -1, 2., 1., colors.ListedColormap(['0.5', 'blue']), 'Ice limit(threshold : emis_LAMB > ' + str(emis_lim_lamb[imo])[0:6] + ')')328 title(str(month[imo]) + ' 2009 - AMSUA 30GHz')329 plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/AMSUA 30_ice_emis_lamb_thresh' + '_' + month[imo] + '2009.png')327 ice_zone_lamb[imo, ilat, ilon] = emis_lamb_month[imo, ilat, ilon] 328 map_cartesian_grid.draw_map_cartes_l(x_coast, y_coast, z_coast, volume, xdist, ydist, ice_zone_lamb[8, :, :], 0.5, 1.02, 0.01, cm.s3pcpn_l_r, 'Sea ice emissivity lamb (threshold : emis_LAMB > ' + str(emis_lim_lamb[imo])[0:6] + ')') 329 title(str(month[imo]) + ' 2009 - AMSUA ' + str(frequ) + 'GHz') 330 plt.savefig('/usr/home/lahlod/twice_d/fig_output_ARCTIC/fig_output_sea_ice_study/ice_class_AMSUA/AMSUA' + str(frequ) + '_ice_emis_lamb_thresh' + '_' + month[imo] + '2009.png') 330 331 331 332 … … 343 344 for imo in range (0, M): 344 345 ice_spec = reshape(ice_zone_spec[imo, :, :], size(ice_zone_spec[imo, :, :])) 345 nb_pix_spec[imo] = len(np.where(i ce_spec == 1.)[0])346 nb_pix_spec[imo] = len(np.where(isnan(ice_spec) == False)[0]) 346 347 total_ice_area_spec[imo] = (pix_area * nb_pix_spec[imo]) + (926. * pix_area) 347 348 … … 353 354 for imo in range (0, M): 354 355 ice_lamb = reshape(ice_zone_lamb[imo, :, :], size(ice_zone_lamb[imo, :, :])) 355 nb_pix_lamb[imo] = len(np.where(i ce_lamb == 1.)[0])356 nb_pix_lamb[imo] = len(np.where(isnan(ice_lamb) == False)[0]) 356 357 total_ice_area_lamb[imo] = (pix_area * nb_pix_lamb[imo]) + (926. * pix_area) 357 358 … … 364 365 ###################### 365 366 print 'start stacking in .dat file' 366 data_classif = open ('/net/argos/data/parvati/lahlod/ARCTIC/AMSU A_ice_class/AMSUA30_data_classification_parameters_ice_no-ice_2009.dat', 'a')367 data_classif = open ('/net/argos/data/parvati/lahlod/ARCTIC/AMSUB_ice_class/AMSUB'+str(frequ)+'_data_classification_parameters_ice_no-ice_2009.dat', 'a') 367 368 for imo in range (0, M): 368 369 for ii in range (0, 50): … … 391 392 for imo in range (0, M): 392 393 print 'stack in file month ' + str(month[imo]) 393 rootgrp = Dataset('/net/argos/data/parvati/lahlod/ARCTIC/AMSU A_ice_class/cartesian_grid_map_ice_no-ice_' + month[imo] + '2009_AMSUA30_spec_lamb_thresholds.nc', 'w', format='NETCDF3_CLASSIC')394 rootgrp = Dataset('/net/argos/data/parvati/lahlod/ARCTIC/AMSUB_ice_class/cartesian_grid_map_ice_no-ice_' + month[imo] + '2009_AMSUB' + str(frequ) + '_spec_lamb_thresholds.nc', 'w', format='NETCDF3_CLASSIC') 394 395 rootgrp.createDimension('longitude', len(xdist)) 395 396 rootgrp.createDimension('latitude', len(ydist))
Note: See TracChangeset
for help on using the changeset viewer.
