;$Id$ @apply_filters @compute_class_cond ;-------------------------------------------------------------------------------------------- ; probabilistic supercooled for 2d naive bayesian with classifier region defined by T_Opa ; ;-------------------------------------------------------------------------------------------- pro prob_sc_2d, topa_btd8511 = topa_btd8511, topa_btd1112 = topa_btd1112, topa_logbt11std = topa_logbt11std, $ etropo_emiss375 = etropo_emiss375, etropo_logcod065 = etropo_logcod065, etropo_ndsi = etropo_ndsi, $ etropo_bt11bt67covari=etropo_bt11bt67covari, etropo_logcod138=etropo_logcod138, $ etropo_btd8511 = etropo_btd8511, etropo_btd1112 = etropo_btd1112, $ etropo_logbt11std = etropo_logbt11std, $ etropo_ref375 = etropo_ref375, etropo_ref160 = etropo_ref160, $ topa_ref375 = topa_ref375, topa_ref160 = topa_ref160, $ ref160_ref138 = ref160_ref138, $ btd1112_btd8511 = btd1112_btd8511, $ day = day, night=night ;------------------------------- ; read in data ;------------------------------- sensor_name='modis' ;sensor_name='viirs' training_source = 'phase_training_modis.sav' restore, training_source f = training_data training_data = !null ;f.surface_elevation = f.surface_elevation * 1000.0 ; km to m ;------------------------------- ; set any geographical limits ;------------------------------- cod_clear_thresh = 0.100000 cod_ice_water_restore_thresh = 0.0 ;0.100000 cod_ice_water_remove_thresh = 0.01 lat_min = -90.0 lat_max = 90.0 etrop_cloud_thresh = 0.9 count_min = 10.0 missing = -999.0 ;----- set solar limits day_string = 'all' solzen_min = 0.0 solzen_max = 180.0 if (keyword_set(day)) then begin day_string = 'day' solzen_min = 0.0 solzen_max = 80.0 endif if (keyword_set(night)) then begin day_string = 'night' solzen_min = 90.0 solzen_max = 180.0 endif nsfc = 7 nc_header = '2D SC LUT' if (keyword_set(etropo_logcod065)) then begin classifier_name = 'etropo_logcod065' xtitle = 'etropo' x = f.emiss_tropo_11_0um_nom x_min = -0.4 & x_bin = 0.025 & nbins_x = 60 ytitle = 'logcod065' idx = where(cod065 gt 0.0,cc) y = f.cod065 y[*] = missing idx = where(cod065 gt 0) & if (cc gt 0) then y[idx] = alog10(cod065[idx]) ;y_min = -20 & y_bin = 1.0 & nbins_y = 100.0 y_min = -2 & y_bin = 0.1 & nbins_y = 40.0 glintzen_min = 40.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 1 endif if (keyword_set(etropo_logcod138)) then begin classifier_name = 'etropo_logcod138' xtitle = 'etropo' x = f.emiss_tropo_11_0um_nom x_min = -0.4 & x_bin = 0.025 & nbins_x = 60 ytitle = 'logcod138' idx = where(cod138 gt 0.0,cc) y = f.cod138 y[*] = missing idx = where(cod138 gt 0) & if (cc gt 0) then y[idx] = alog10(cod138[idx]) y_min = -2 & y_bin = 0.1 & nbins_y = 40.0 glintzen_min = 0.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(etropo_emiss375)) then begin classifier_name = 'etropo_emiss375' xtitle = 'etropo' x = f.emiss_tropo_11_0um_nom x_min = -0.4 & x_bin = 0.025 & nbins_x = 60 ytitle = 'emiss375' y = emiss375 y_min = 0.25 & y_bin = 0.05 & nbins_y = 70.0 glintzen_min = 40.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(topa_ref375)) then begin classifier_name = 'topa_ref375' xtitle = 'Opaque Cloud Temperature' x = f.cld_temp_opaque x_min = 180.0 & x_bin = 5.0 & nbins_x = 32 ytitle = 'ref375' y = f.ref_375um y_min = -20.0 & y_bin = 0.5 & nbins_y = 150.0 glintzen_min = 40.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(etropo_ref375)) then begin classifier_name = 'etropo_ref375' xtitle = 'etropo' x = f.emiss_tropo_11_0um_nom x_min = -0.4 & x_bin = 0.025 & nbins_x = 60 ytitle = 'ref375' y = f.ref_375um y_min = -20.0 & y_bin = 0.5 & nbins_y = 50.0 glintzen_min = 40.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(etropo_ref160)) then begin classifier_name = 'etropo_ref160' xtitle = 'etropo' x = f.emiss_tropo_11_0um_nom x_min = -0.4 & x_bin = 0.025 & nbins_x = 60 ytitle = 'ref160' y = f.ref_160um y_min = -2.0 & y_bin = 1.0 & nbins_y = 80.0 glintzen_min = 40.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(topa_ref160)) then begin classifier_name = 'topa_ref160' xtitle = 'Opaque Cloud Temperature' x = f.cld_temp_opaque x_min = 180.0 & x_bin = 5.0 & nbins_x = 32 ytitle = 'ref160' y = f.ref_160um y_min = -2.0 & y_bin = 1.0 & nbins_y = 80.0 glintzen_min = 40.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(etropo_bt11bt67covari)) then begin classifier_name = 'etropo_bt11bt67covari' xtitle = 'etropo' x = f.emiss_tropo_11_0um_nom x_min = -0.4 & x_bin = 0.025 & nbins_x = 60 ytitle = 'bt11bt67covar' y = f.bt11bt67covar y_min = -2 & y_bin = 0.2 & nbins_y = 50.0 glintzen_min = 0.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(etropo_ndsi)) then begin classifier_name = 'etropo_ndsi' xtitle = 'etropo' x = f.emiss_tropo_11_0um_nom x_min = -0.4 & x_bin = 0.025 & nbins_x = 60 ytitle = 'ndsi' y = f.emiss_375um y_min = 0.25 & y_bin = 0.05 & nbins_y = 70.0 glintzen_min = 40.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(etropo_btd8511)) then begin classifier_name = 'etropo_btd8511' xtitle = 'etropo' x = f.emiss_tropo_11_0um_nom x_min = -0.4 & x_bin = 0.025 & nbins_x = 60 ytitle = 'Btd 8.5 - 11' x = f.temp_8_5um_nom - f.temp_11_0um_nom y_min = -5 & y_bin = 0.25 & nbins_y = 60.0 glintzen_min = 0.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(btd1112_btd8511)) then begin classifier_name = 'btd1112_btd8511' xtitle = 'Btd 11 - 12' x = f.temp_11_0um_nom - f.bt_12um x_min = -2 & x_bin = 0.25 & nbins_x = 40.0 ytitle = 'Btd 8.5 - 11' y = f.temp_8_5um_nom - f.temp_11_0um_nom y_min = -5 & y_bin = 0.25 & nbins_y = 60.0 glintzen_min = 0.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(topa_btd8511)) then begin classifier_name = 'topa_btd8511' xtitle = 'Opaque Cloud Temperature' x = f.cld_temp_opaque x_min = 180.0 & x_bin = 5.0 & nbins_x = 32 ytitle = 'Btd 8.5 - 11' y = f.temp_8_5um_nom - f.temp_11_0um_nom y_min = -5 & y_bin = 0.25 & nbins_y = 60.0 glintzen_min = 0.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(etropo_btd1112)) then begin classifier_name = 'etropo_btd1112' xtitle = 'etropo' x = f.emiss_tropo_11_0um_nom x_min = -0.4 & x_bin = 0.025 & nbins_x = 60 ytitle = 'Btd 11 - 12' y = f.temp_11_0um_nom - f.bt_12um y_min = -2 & y_bin = 0.25 & nbins_y = 40.0 glintzen_min = 0.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(topa_btd1112)) then begin classifier_name = 'topa_btd1112' xtitle = 'Opaque Cloud Temperature' x = f.cld_temp_opaque x_min = 180.0 & x_bin = 5.0 & nbins_x = 32 ytitle = 'Btd 11 - 12' y = f.temp_11_0um_nom - f.bt_12um y_min = -2 & y_bin = 0.25 & nbins_y = 40.0 glintzen_min = 0.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(etropo_logbt11std)) then begin classifier_name ='etropo_logbt11std' xtitle = 'etropo' x = f.emiss_tropo_11_0um_nom x_min = -0.4 & x_bin = 0.025 & nbins_x = 60 x_min = -0.6 & x_bin = 0.025 & nbins_x = 80 ytitle = 'log10 BT11 Stddev' y = logt11_std y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0 glintzen_min = 0.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(topa_logbt11std)) then begin classifier_name ='topa_logbt11std' xtitle = 'Opaque Cloud Temperature' x = f.cld_temp_opaque x_min = 180.0 & x_bin = 5.0 & nbins_x = 32 ytitle = 'log10 BT11 Stddev' y = f.logt11_std y_min = -1.5 & y_bin = 0.25 & nbins_y = 15.0 glintzen_min = 0.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif if (keyword_set(ref160_ref138)) then begin classifier_name = 'ref160_ref138' xtitle = 'ref160' x = f.ref_160um x_min = -2.0 & x_bin = 1.0 & nbins_x = 65 ytitle = 'ref138' y = f.ref_138um y_min = -2.0 & y_bin = 0.5 & nbins_y = 80.0 glintzen_min = 0.0 glintzen_max = 180.0 zsfc_min = -100.0 zsfc_max = 20000.0 snow_class_min = 1 snow_class_max = 3 endif ;------------------------------------------------ ; make bin bounds ;------------------------------------------------ x_x = x_min + findgen(nbins_x)*x_bin + x_bin/2 x_y = y_min + findgen(nbins_y)*y_bin + y_bin/2 ;----------------------------------------------- ; make surface type ;---------------------------------------------- ;sst_back_uni = make_array(n_elements(lat), /FLOAT, VALUE=0.0) ;bayes_sfc_type_1d, land_class, snow_class, sst_back_uni, lon, lat, EMISS375_SURFACE, sfc_type_mask sfc_type_mask = f.bayes_mask_sfc_type ;----------------------------------------------- ; apply filters ;---------------------------------------------- apply_filters, f.CLOSEST_CALIPSO_CLOUD_PHASE_TOP_LAYER, f.CLOSEST_CALIPSO_CLOUD_FRACTION_5KM, $ f.CLOSEST_CALIPSO_COD, f.CLOSEST_CALIPSO_COD_ICE, $ cod_clear_thresh, cod_ice_water_restore_thresh, cod_ice_water_remove_thresh, $ f.solar_zenith_angle, solzen_min, solzen_max, $ f.land_class, f.glint_zenith_angle, glintzen_min, glintzen_max, $ f.snow_class, snow_class_min, snow_class_max, $ f.surface_elevation, zsfc_min, zsfc_max, $ f.latitude, lat_min, lat_max, $ f.emiss_tropo_11_0um_nom, etrop_cloud_thresh, missing ;------------------------------- ; select valid data within these limits ;------------------------------- idx = where(f.CLOSEST_CALIPSO_CLOUD_PHASE_TOP_LAYER ge 0.0 and x ne missing and y ne missing,cc) print, 'total number after filtering= ', cc ;--- through away filtered data x = x[idx] y = y[idx] tc_lidar = f[idx].closest_calipso_top_mid_temperature tc_opa = f[idx].cld_temp_opaque zc_opa = f[idx].cld_height_opaque cod_lidar = f[idx].CLOSEST_CALIPSO_COD etropo = f[idx].EMISS_TROPO_11_0UM_NOM phase_lidar = f[idx].CLOSEST_CALIPSO_CLOUD_PHASE_TOP_LAYER sfc_type_mask = sfc_type_mask[idx] btd8511 = f[idx].temp_8_5um_nom - f[idx].temp_11_0um_nom ;---------------------------------------------------------- ; change sc_mask_true = 0 = no, 1 = yes ;---------------------------------------------------------- sc_mask_true = phase_lidar ;---- prob relative to all pixels ;sc_mask_true[*] = -1 ;idx = where(phase_lidar eq 0 or phase_lidar eq 1 or phase_lidar eq 3 or tc_lidar gt 273.15, cc) ;if (cc gt 0) then sc_mask_true[idx] = 0 ;idx = where(phase_lidar eq 2 and tc_lidar gt 0.0 and tc_lidar le 273.15, cc) ;if (cc gt 0) then sc_mask_true[idx] = 1 ;---- prob relative to all clouds sc_mask_true[*] = -1 idx = where(phase_lidar eq 1 or phase_lidar eq 3 or tc_lidar gt 273.15 or tc_lidar lt 233.0, cc) if (cc gt 0) then sc_mask_true[idx] = 0 idx = where(phase_lidar eq 2 and tc_lidar gt 0.0 and tc_lidar le 273.15, cc) if (cc gt 0) then sc_mask_true[idx] = 1 ;------------------------------------------------------- ; add limits to exclude pixels from analysis ;------------------------------------------------------- idx = where(tc_opa gt 273.15,cc) if (cc gt 0) then sc_mask_true[idx] = -1 idx = where(tc_opa lt 233.0,cc) if (cc gt 0) then sc_mask_true[idx] = -1 ;idx = where(cod_lidar lt 1.0,cc) ;if (cc gt 0) then sc_mask_true[idx] = -1 idx = where(btd8511 gt 2.0,cc) if (cc gt 0) then sc_mask_true[idx] = -1 idx = where(etropo lt 0.1,cc) if (cc gt 0) then sc_mask_true[idx] = -1 ;---- prob relative to water clouds ;sc_mask_true[*] = -1 ;idx = where(phase_lidar eq 1 or phase_lidar eq 3, cc) ;if (cc gt 0) then sc_mask_true[idx] = -1 ;idx = where(phase_lidar eq 2 and tc_lidar gt 273.15, cc) ;if (cc gt 0) then sc_mask_true[idx] = 0 ;idx = where(phase_lidar eq 2 and tc_lidar gt 0.0 and tc_lidar le 273.15, cc) ;if (cc gt 0) then sc_mask_true[idx] = 1 laplace_flag = 0 nb_rank = 2 ;------------------------------------------------------------- ; make arrays ;------------------------------------------------------------- cond_ratio_table = make_array(nbins_x, nbins_y, nsfc, /FLOAT, VALUE=MISSING) obs_prob_table = make_array(nbins_x, nbins_y, nsfc, /FLOAT, VALUE=0.0) sc_fraction = make_array(nsfc,/FLOAT,VALUE=missing) observation_count = make_array(nsfc,/FLOAT,VALUE=missing) ;-------------------------------------------------------------------- ; loop over surface type ;-------------------------------------------------------------------- for isfc = 0, nsfc - 1 do begin idx = where(sfc_type_mask eq isfc + 1 and sc_mask_true ne -1 ,cc) if (cc gt 0) then begin x_sfc = x[idx] y_sfc = y[idx] sc_mask_true_sfc = sc_mask_true[idx] observation_count[isfc] = cc print, 'isfc, count = ', isfc, cc endif else begin goto, skip_sfc endelse ;--- make priors idx = where(sc_mask_true_sfc ge 0, count_all) idx = where(sc_mask_true_sfc eq 1, count_sc) sc_fraction[isfc] = float(count_sc) / float(count_all) prior_yes = sc_fraction[isfc] prior_no = 1.0 - prior_yes print, 'prior yes = ', prior_yes cond_ratio = make_array(nbins_x, nbins_y, /FLOAT, VALUE=MISSING) post_prob = make_array(nbins_x, nbins_y, /FLOAT, VALUE=MISSING) obs_prob = make_array(nbins_x, nbins_y, /FLOAT, VALUE=0.0) idx = where(sc_mask_true_sfc ge 0 and sc_mask_true_sfc eq 1,count_yes_all) idx = where(sc_mask_true_sfc ge 0 and sc_mask_true_sfc eq 0,count_no_all) cond_ratio[*] = 1.0 max_r = 100.0 for i = 0, nbins_x -1 do begin x_min_bin = x_min + i*x_bin x_max_bin = x_min + (i+1)*x_bin for j = 0, nbins_y -1 do begin y_min_bin = y_min + j*y_bin y_max_bin = y_min + (j+1)*y_bin idx_all = where(x_sfc ge x_min_bin and x_sfc le x_max_bin and $ y_sfc ge y_min_bin and y_sfc le y_max_bin and $ sc_mask_true_sfc ne -1, count_all_bin) ; print, x_min_bin, x_max_bin, y_min_bin, y_max_bin ;if (count_all_bin gt 1) then stop idx_yes = where(sc_mask_true_sfc[idx_all] eq 1, count_yes_bin) idx_no = where(sc_mask_true_sfc[idx_all] eq 0, count_no_bin) compute_class_cond, nb_rank, laplace_flag, $ count_min, count_yes_all, count_no_all, count_yes_bin, count_no_bin, $ max_r, prior_yes, cond_yes, cond_no, r, posterior_prob cond_ratio[i,j] = r post_prob[i,j] = posterior_prob obs_prob[i,j] = float(count_all_bin) ; print, i, j, count_all_bin, count_yes_bin, count_no_bin endfor endfor ;--- normalize obs_prob idx = where(obs_prob gt 0,cc) if (cc gt 0) then obs_prob[idx] = obs_prob[idx] / total(obs_prob) ;----------------- smooth for ix = 1, nbins_x - 2 do begin for iy = 1, nbins_y - 2 do begin v = obs_prob[ix-1:ix+1,iy-1:iy+1] idx = where(v gt 0,cc) if (cc gt 0) then begin v = cond_ratio[ix-1:ix+1,iy-1:iy+1] cond_ratio[ix,iy] = mean(v[idx]) v = post_prob[ix-1:ix+1,iy-1:iy+1] post_prob[ix,iy] = mean(v[idx]) endif endfor endfor ;------------------- ; fill in holes ;------------------- idx = where(obs_prob eq 0,cc) & print, 'number with zero = ', cc for i = 0, 5 do begin for ix = 1, nbins_x - 2 do begin for iy = 1, nbins_y - 2 do begin if (obs_prob[ix,iy] eq 0) then begin v = obs_prob[ix-1:ix+1,iy-1:iy+1] idx = where(v gt 0,cc) if (cc gt 0) then begin v = cond_ratio[ix-1:ix+1,iy-1:iy+1] cond_ratio[ix,iy] = mean(v[idx]) v = post_prob[ix-1:ix+1,iy-1:iy+1] post_prob[ix,iy] = mean(v[idx]) endif endif endfor endfor endfor ;------------------------------------------------------------------------------------ ; apply table ;------------------------------------------------------------------------------------ n = n_elements(sc_mask_true_sfc) sc_retrieved = make_array(n, /FLOAT, VALUE = missing) for k = 0, n - 1 do begin i = min([nbins_x -1,max([0,fix((x_sfc[k] - x_min) / x_bin)])]) j = min([nbins_y -1,max([0,fix((y_sfc[k] - y_min) / y_bin)])]) if(post_prob[i,j] ne missing) then sc_retrieved[k] = 0 if(post_prob[i,j] ge 0.50) then sc_retrieved[k] = 1 endfor idx = where(sc_mask_true_sfc ge 0 and sc_retrieved ge 0, count_all) idx = where(sc_mask_true_sfc eq 0 and sc_retrieved eq 0, count_agree_nosc) idx = where(sc_mask_true_sfc eq 0 and sc_retrieved eq 1, count_false_sc) idx = where(sc_mask_true_sfc eq 1 and sc_retrieved eq 0, count_missed_sc) idx = where(sc_mask_true_sfc eq 1 and sc_retrieved eq 1, count_agree_sc) idx = where(sc_mask_true_sfc eq 1, count_sc_true) idx = where(sc_retrieved eq 1, count_sc_new) pod_sc = float(count_agree_sc + count_agree_nosc) / count_all missed_sc = float(count_missed_sc) / count_all false_sc = float(count_false_sc) / count_all frac_sc_true = float(count_sc_true) / count_all frac_sc_new = float(count_sc_new) / count_all print, 'frac sc true = ', frac_sc_true print, 'frac sc new = ', frac_sc_new print, 'pod sc = ', pod_sc print, 'false sc = ', false_sc print, 'missed sc = ', missed_sc ;------------------------------------------------------------------------------------ ; make images ;------------------------------------------------------------------------------------ window, 0, xsize = 600, ysize = 600 load_color_table, 33, 0, 1.0, 249 erase, color = 0 view, post_prob, coltab=33, /current, position = [0.10,0.06,0.98,0.98], cbar_position=cbar_position, /nocbar, missing_color=1, $ sds_min=-101, sds_max=1, /mask_missing, title = 'SUPERCOOLED' plot, [0.0,0.0],[100.0,100.0],color = 3, background = 1 , $ xrange = [min(x_x),max(x_x)],ystyle = 1, /nodata, /noerase, $ yrange = [min(x_y),max(x_y)],xstyle = 1, $ xtitle = xtitle, ytitle = ytitle legend_local, ['SUPERCOOLED'], color = 3, /top, /right, charsize = 2, /clear legend_local, [string(isfc+1,format="(I02)")], color = 3, /top, /left, charsize = 2, /clear saveimage, classifier_name+'_sc_prob_lidar_sfc'+string(isfc,format="(I01)")+'.png',/png obs_prob_table[*,*,isfc] = obs_prob cond_ratio_table[*,*,isfc] = cond_ratio skip_sfc: endfor ;sfc type loop ;---------------------------------------------- ; write to netcdf file ;---------------------------------------------- nc_filename = sensor_name + "_" + day_string+"_"+classifier_name+"_sc_prob.nc" sd_id = ncdf_create(nc_filename,/CLOBBER) ncdf_attput,sd_id, "description", nc_header, /CHAR, /global ncdf_attput,sd_id, "sensor", sensor_name, /CHAR, /global ncdf_attput,sd_id, "training_source", training_source, /CHAR, /global ncdf_attput,sd_id, "timestamp",timestamp(), /CHAR, /global ncdf_attput,sd_id, "rank", nb_rank, /SHORT, /global ncdf_attput,sd_id, "x_name", xtitle, /CHAR, /global ncdf_attput,sd_id, "nbins_x", float(nbins_x), /FLOAT, /global ncdf_attput,sd_id, "x_min", float(x_min), /FLOAT, /global ncdf_attput,sd_id, "x_bin", float(x_bin), /FLOAT, /global ncdf_attput,sd_id, "y_name", ytitle, /CHAR, /global ncdf_attput,sd_id, "nbins_y", float(nbins_y), /FLOAT, /global ncdf_attput,sd_id, "y_min", float(y_min), /FLOAT, /global ncdf_attput,sd_id, "y_bin", float(y_bin), /FLOAT, /global ncdf_attput,sd_id, "nsfc", float(nsfc), /FLOAT, /global ncdf_attput,sd_id, "latitude_minimum", lat_min, /FLOAT, /global ncdf_attput,sd_id, "latitude_maximum", lat_max, /FLOAT, /global ncdf_attput,sd_id, "calipso_cod_clear_threshold", cod_clear_thresh, /FLOAT, /global ncdf_attput,sd_id, "solar_zenith_minimum", solzen_min, /FLOAT, /global ncdf_attput,sd_id, "solar_zenith_maximum", solzen_max, /FLOAT, /global ncdf_attput,sd_id, "glint_zenith_minimum", glintzen_min, /FLOAT, /global ncdf_attput,sd_id, "glint_zenith_maximum", glintzen_max, /FLOAT, /global ncdf_attput,sd_id, "surface_elevation_minimum", zsfc_min, /FLOAT, /global ncdf_attput,sd_id, "surface_elevation_maximum", zsfc_max, /FLOAT, /global ncdf_attput,sd_id, "snow_class_minimum", snow_class_min, /FLOAT, /global ncdf_attput,sd_id, "snow_class_maximum", snow_class_max, /FLOAT, /global ;--- get dimensions nx_id = NCDF_DIMDEF(sd_id, 'nbins_x', nbins_x) ny_id = NCDF_DIMDEF(sd_id, 'nbins_y', nbins_y) nk_id = NCDF_DIMDEF(sd_id, 'nsfc', nsfc) ;--- create variable sds_id1 = NCDF_VARDEF(sd_id, 'cond_ratio_table', [nx_id,ny_id,nk_id], /FLOAT) sds_id2 = NCDF_VARDEF(sd_id, 'obs_prob_table', [nx_id,ny_id,nk_id], /FLOAT) sds_id3 = NCDF_VARDEF(sd_id, 'observation_count', [nk_id], /FLOAT) sds_id4 = NCDF_VARDEF(sd_id, 'prior_yes', [nk_id], /FLOAT) NCDF_CONTROL, sd_id, /ENDEF ; Put the file into data mode NCDF_VARPUT, sd_id, sds_id1, cond_ratio_table NCDF_VARPUT, sd_id, sds_id2, obs_prob_table NCDF_VARPUT, sd_id, sds_id3, observation_count NCDF_VARPUT, sd_id, sds_id4, sc_fraction ;--- close netcdf files ncdf_close, sd_id end