;$Id$
@apply_filters
@fmft
@write_2d_table
@make_true_phase
@smooth_nan
@ecm_sfc_type
;--------------------------------------------------------------------------------------------
;  probabilistic phase for 2d naive bayesian
;
;--------------------------------------------------------------------------------------------
pro prob_mask_phase_2d, day = day, night=night, sensor_name = sensor_name, noplot=noplot, new=new, old=old, $
                        test_name=test_name, sfc_type_on_flag_input = sfc_type_on_flag_input, $
                        path_dir = path_dir, nomod = nomod, Time_Diff_Max = Time_Diff_Max, $
                        cod_clear_ice_thresh = cod_clear_ice_thresh, cod_clear_water_thresh = cod_clear_water_thresh, yesplot=yesplot, $
			isccp = isccp

if (~keyword_set(path_dir)) then path_dir='./'
if (~keyword_set(sensor_name)) then sensor_name='myd02ssh'
if (~keyword_set(Time_Diff_Max)) then Time_Diff_Max = 10.0  ; minutes
if (~keyword_set(cod_clear_ice_thresh)) then cod_clear_ice_thresh = 0.1
if (~keyword_set(cod_clear_water_thresh)) then cod_clear_water_thresh = 0.1
if (sensor_name eq 'myd02ssh') then training_source = 'ecm_training_myd02ssh_train.sav'
if (sensor_name eq 'viirs') then training_source = 'ecm_training_viirs_train_snpp_5+333_2013.sav'
if (sensor_name eq 'abhi') then training_source = 'ecm_training_abhi_train.sav'
if (sensor_name eq 'avhrr') then training_source = 'ecm_training_avhrr_train.sav'
if (sensor_name eq 'mtsat2') then training_source = 'ecm_training_mtsat2_train.sav'
if (sensor_name eq 'vgac') then training_source = 'ecm_training_vgac_train.sav'
if (sensor_name eq 'seviri') then training_source = 'ecm_training_seviri_train.sav'

;--- check for sensor and test inconsistencies
if (strpos(sensor_name,'viirs') lt 0 and strpos(test_name,'dnb') gt 0) then begin
   print, 'DNB classifiers impossible'
   return
endif

;-- control training data modification in make_true_phase routine
mod_flag = 1
if (keyword_set(nomod)) then mod_flag = 0

;-------------------------------
; read in data
;-------------------------------
training_path = './training_data/'
restore, training_path + training_source
f = training_data
training_data = !null
ntrain = f.length

nb_rank = 2

nsfc = 7

if (~keyword_set(sfc_type_on_flag_input)) then sfc_type_on_flag_input = make_array(nsfc,/INT,VALUE=1)

;---- some parameters
laplace_flag = 0
max_r = 1000.0
count_min = 10
missing = -999.0
rut_solzen_thresh = 80.0

;--- partial cloud filter
cfrac_min = 0.1
cfrac_max = 0.9

;--- set any geographical limits
lat_min = -90.0
lat_max = 90.0

;---- adjust nasa dnb to look like noaa dnb [REMOVE WHEN ADJUSTMENT in TRAINING DATA]
Dnb_Coef = [-0.373685,0.977945,-0.000261637]
Dnb_Coef = [-0.118767,0.962452,-0.000144502] ; Denis

idx = where(f.refl_lunar_dnb_nom ne missing,cc)
if (cc gt 0) then begin
   f[idx].refl_lunar_dnb_nom = Dnb_Coef[0] + $
                              Dnb_Coef[1]*f[idx].refl_lunar_dnb_nom + $
                              Dnb_Coef[2]*f[idx].refl_lunar_dnb_nom^2
endif

;--- product filtering
cod_cloud_thresh = 1.0
etrop_cloud_thresh = 0.9

;---  default classifier limits - only modify in classifier definition if needed
zen_min = 0.0
zen_max = 90.0
solzen_min = 0.0
solzen_max = 180.0
solzen_day_max = 85.0
solglintzen_min = 0.0
solglintzen_max = 180.0
solscatang_min = 0.0
solscatang_max = 180.0
solglint_mask_min = 0
solglint_mask_max = 1
lunzen_min = 0.0
lunzen_max = 180.0
lunglintzen_min = 0.0
lunglintzen_max = 180.0
lunscatang_min = 0.0
lunscatang_max = 180.0
lunglint_mask_min = 0
lunglint_mask_max = 1
coast_mask_min = 0
coast_mask_max = 1
city_mask_min = 0
city_mask_max = 1
moon_illum_frac_min = 0.0
moon_illum_frac_max = 100.0
slant_tpw_min = 0.0
slant_tpw_max = 100.0
tsfc_min = 180.0
tsfc_max = 340.0
zsfc_min = -100.0
zsfc_max = 20000.0
zsfc_std_min = 00.0
zsfc_std_max = 1000.0
snow_class_min = 1
snow_class_max = 3
sfc_type_on_flag = sfc_type_on_flag_input
zsfc_std_max_iruni = 300.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 = solzen_day_max
endif
if (keyword_set(night)) then begin
   day_string = 'night'
   solzen_min = 90.0
   solzen_max = 180.0
endif

;------------------------------------------------------------------------------
; small positive (near zero) dnb reflectances should be missing
;------------------------------------------------------------------------------
idx = where(f.refl_lunar_dnb_nom lt 0.1,cc) & if (cc gt 0) then f[idx].refl_lunar_dnb_nom = missing

;------------------------------------------------------------------------------
; loop through classifiers
;------------------------------------------------------------------------------
nc_header = '2D Mask and Phase LUT'

print, '==============================================='
print , "2D Classifier =  ", test_name

case test_name of

"refldnb_refl38": begin
   classifier_name = 'refldnb_refl38'
   xtitle = 'refldnb'
   x = f.refl_lunar_dnb_nom
   x_min = -5 & x_bin = 1.0 & nbins_x = 80
   ytitle = 'refl38'
   y = f.refl_3_75um_nom  & y_min = -20.0 & y_bin = 1.0 & nbins_y = 50.0
   nchan_used = 2 & channel_wvl_used = [650,3750]
   tsfc_min = 240.0
   solzen_min = 90.0
   solzen_max = 180.0
   lunzen_min = 0.0
   lunzen_max = 80.0
   lunglint_mask_min = 0
   lunglint_mask_max = 0
   moon_illum_frac_min = 25.0
   moon_illum_frac_max = 100.0
end

"refldnb_btd3811": begin
   classifier_name = 'refldnb_btd3811'
   xtitle = 'refldnb'
   x = f.refl_lunar_dnb_nom
   x_min = -5 & x_bin = 1.0 & nbins_x = 80
   ytitle = 'btd3811'
   y = f.temp_3_75um_nom - f.temp_11_0um_nom
   y_min = -10 & y_bin = 1.0 & nbins_y = 40.0
   nchan_used = 2 & channel_wvl_used = [650,3750]
   tsfc_min = 240.0
   solzen_min = 90.0
   solzen_max = 180.0
   lunzen_min = 0.0
   lunzen_max = 80.0
   lunglint_mask_min = 0
   lunglint_mask_max = 0
   moon_illum_frac_min = 25.0
   moon_illum_frac_max = 100.0
end

"logzopa_refldnb": begin
   classifier_name = 'logzopa_refldnb'
   xtitle = 'logzopa' & x = f.cld_height_opaque & x_min = 1.5 & x_bin = 0.1 & nbins_x = 30.0
   idx = where(x gt 0.0,cc) & if (cc gt 0) then x[idx] = alog10(x[idx])
   idx = where(x eq 0.0,cc) & if (cc gt 0) then x[idx] = x_min
   ytitle = 'refldnb' & y = f.refl_lunar_dnb_nom & y_min = -5 & y_bin = 2.5 & nbins_y = 50
   nchan_used = 2 & channel_wvl_used = [650,1100]
   tsfc_min = 240.0
   solzen_min = 90.0
   solzen_max = 180.0
   lunzen_min = 0.0
   lunzen_max = 80.0
   lunglint_mask_min = 0
   lunglint_mask_max = 0
   moon_illum_frac_min = 25.0
   moon_illum_frac_max = 100.0
end

"etropo11_logzopa": begin
   classifier_name = 'etropo11_logzopa'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom &  x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'logzopa' & y = f.cld_height_opaque & y_min = 0.5 & y_bin = 0.15 & nbins_y = 30.0
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [11000]
end

"etropo11_logcod065": begin
   classifier_name = 'etropo11_logcod065'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom &  x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'logcod065'
   y = f.cld_opd_mask_0_65um_nom
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   y_min = -2 & y_bin = 0.1 & nbins_y = 40.0
   nchan_used = 2 & channel_wvl_used = [650,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solglint_mask_min = 0
   solglint_mask_max = 0
   snow_class_min = 1
   snow_class_max = 1
end

"etropo11_logcod138": begin
   classifier_name = 'etropo11_logcod138'
   xtitle = 'etropo11'
   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)
   yy = f.cld_opd_mask_1_38um_nom
   y[*] = missing
   idx = where(yy gt 0) & if (cc gt 0) then y[idx] = alog10(yy[idx])
   y_min = -2 & y_bin = 0.1 & nbins_y = 40.0
   yy = !null
   nchan_used = 2 & channel_wvl_used = [1380,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   slant_tpw_min = 0.5
end

"etropo11_emiss3811_day": begin
   classifier_name = 'etropo11_emiss3811_day'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'emiss3811'
   y = f.emiss_3_75um_nom
   y_min = 0.25 & y_bin = 0.05 & nbins_y = 70.0
   nchan_used = 2 & channel_wvl_used = [3750,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solscatang_min = 95.0
   solglint_mask_min = 0
   solglint_mask_max = 0
   tsfc_min = 240.0
end

"etropo11_emiss3811_night": begin
   classifier_name = 'etropo11_emiss3811_night'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'emiss3811'
   y = f.emiss_3_75um_nom
   y_min = 0.25 & y_bin = 0.05 & nbins_y = 70.0
   nchan_used = 2 & channel_wvl_used = [3750,11000]
   solzen_min = 90.0
   solzen_max = 180.0
   solglint_mask_min = 0
   solglint_mask_max = 0
   tsfc_min = 240.0
end

"etropo10_emiss3810": begin
   classifier_name = 'etropo10_emiss3810'
   xtitle = 'etropo10'
   x = f.emiss_tropo_10_4um_nom
   x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'emiss3810'
   y = emiss3810
   y = f.emiss_3_75um_nom
   y_min = 0.25 & y_bin = 0.05 & nbins_y = 70.0
   nchan_used = 2 & channel_wvl_used = [3750,10400]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solglint_mask_min = 0
   solglint_mask_max = 0
   tsfc_min = 240.0
end

"etropo11_refl38": begin
   classifier_name = 'etropo11_refl38'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom & x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'refl38'
   y = f.refl_3_75um_nom  & y_min = -20.0 & y_bin = 0.5 & nbins_y = 50.0
   nchan_used = 2 & channel_wvl_used = [3750,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solglint_mask_min = 0
   solglint_mask_max = 0
   tsfc_min = 240.0
end

"etropo11_refl160": begin
   classifier_name = 'etropo11_refl160'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'refl160'
   y = f.refl_1_60um_nom
   y_min = -2.0 & y_bin = 1.0 & nbins_y = 80.0
   nchan_used = 2 & channel_wvl_used = [1600,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solglint_mask_min = 0
   solglint_mask_max = 0
end

"bt11_refl065": begin
   classifier_name = 'bt11_refl065'
   xtitle = 'bt11' & x = f.temp_11_0um_nom & x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'refl065' & y = f.refl_0_65um_nom & y_min = -10.0 & y_bin = 2.0 & nbins_y = 50
   nchan_used = 2 & channel_wvl_used = [650,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solglint_mask_min = 0
   solglint_mask_max = 0
   snow_class_max = 1
end

"etropo11_bt6711covar": begin
   classifier_name = 'etropo11_bt6711covar'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'bt6711covar'
   y = f.temp_11um_vs_67um_covar_5x5
   y_min = -4 & y_bin = 0.15 & nbins_y = 80
   nchan_used = 2 & channel_wvl_used = [6700,11000]
end

"zopa_bt6711covar": begin
   classifier_name = 'zopa_bt6711covar'
   xtitle = 'zopa'
   x = f.cld_height_opaque
   x_min = 0.0 & x_bin = 1000.0 & nbins_x = 20
   ytitle = 'bt6711covar'
   y = f.temp_11um_vs_67um_covar_5x5
   y_min = -2 & y_bin = 0.2 & nbins_y = 50.0
   nchan_used = 2 & channel_wvl_used = [6700,11000]
end

"logzopa_btd1112": begin
   classifier_name = 'logzopa_btd1112'
   xtitle = 'logzopa' & x = f.cld_height_opaque & x_min = 1.5 & x_bin = 0.1 & nbins_x = 30.0
   idx = where(x gt 0.0,cc) & if (cc gt 0) then x[idx] = alog10(x[idx])
   idx = where(x eq 0.0,cc) & if (cc gt 0) then x[idx] = x_min
   ytitle = 'btd1112' & y = f.temp_11_0um_nom - f.temp_12_0um_nom & y_min = -2.0 & y_bin = 0.2 & nbins_y = 50
   nchan_used = 2 & channel_wvl_used = [11000,12000]
   coast_mask_max = 0
;  zsfc_max = 2000.0
   tsfc_min = 240.0
end

"dzopasfc_logcod065": begin
   classifier_name = 'dzopasfc_logcod065'
   xtitle = 'dzopasfc'
   x = f.cld_height_opaque - f.surface_elevation &  x_min = -2000.0 & x_bin = 200.0 & nbins_x = 90
   ytitle = 'logcod065'
   y = f.cld_opd_mask_0_65um_nom
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   y_min = -2 & y_bin = 0.1 & nbins_y = 40.0
   nchan_used = 2 & channel_wvl_used = [650,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solglint_mask_min = 0
   solglint_mask_max = 0
   snow_class_min = 1
   snow_class_max = 1
end

"dzopasfc_drefl065clr": begin
   classifier_name = 'dzopasfc_logcod065'

   xtitle = 'dzopasfc'
   x = f.cld_height_opaque - f.surface_elevation &  x_min = -2000.0 & x_bin = 200.0 & nbins_x = 90

   ytitle = 'drefl065clr'
   y = f.refl_0_65um_nom - f.refl_0_65um_nom_clear_sky & y_min = -40.0 & y_bin = 1.0 & nbins_y = 120

   nchan_used = 2 & channel_wvl_used = [650,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solglint_mask_min = 0
   solglint_mask_max = 0
   snow_class_min = 1
   snow_class_max = 1
end

"dzopasfc_ndsi": begin
   classifier_name = 'dzopasfc_ndsi'

   xtitle = 'dzopasfc'
   x = f.cld_height_opaque - f.surface_elevation &  x_min = -2000.0 & x_bin = 200.0 & nbins_x = 90

   ytitle = 'ndsi'
   y = ndsi(f.refl_0_65um_nom, f.refl_1_60um_nom, missing) & y_min = -0.5 & y_bin = 0.025 & nbins_y = 80

   nchan_used = 2 & channel_wvl_used = [650,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solglint_mask_min = 0
   solglint_mask_max = 0
   snow_class_min = 1
   snow_class_max = 3
end


"etropo11_ndsi": begin
   classifier_name = 'etropo11_ndsi'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.4 & x_bin = 0.025 & nbins_x = 60

;  ytitle = 'ndsi'
;  y = (f.refl_0_65um_nom - f.refl_1_60um_nom) / $
;      (f.refl_0_65um_nom + f.refl_1_60um_nom)
;  idx = where(f.refl_0_65um_nom eq missing or f.refl_1_60um_nom eq missing,cc)
;  if (cc gt 0) then y[idx] = missing
;  y_min = -0.5 & y_bin = 0.025 & nbins_y = 60.0

   ytitle = 'ndsi'
   y = ndsi(f.refl_0_65um_nom, f.refl_1_60um_nom, missing) & y_min = -0.5 & y_bin = 0.025 & nbins_y = 60

   nchan_used = 3 & channel_wvl_used = [650,1600,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solglint_mask_min = 0
   solglint_mask_max = 0
   coast_mask_max = 0
end
 
"bt11ratio_btd8511": begin
   classifier_name = 'bt11ratio_btd8511'
   xtitle = 'bt11ratio'
   x = btratio(f.temp_11_0um_nom,f.surface_temperature_nwp,f.tropopause_temperature_nwp,missing)
   x_min = -0.1 & x_bin = 0.02 & nbins_x = 60
   ytitle = 'btd8511'
   y = f.temp_8_5um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.25 & nbins_y = 60.0
   nchan_used = 2 & channel_wvl_used = [8500,11000]
end

"etropo11_btd8511": begin
   classifier_name = 'etropo11_btd8511'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'btd8511'
   y = f.temp_8_5um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.25 & nbins_y = 60.0
   nchan_used = 2 & channel_wvl_used = [8500,11000]
   tsfc_min = 230.0
end

"etropo11_btd8511_night": begin
   classifier_name = 'etropo11_btd8511_night'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'btd8511'
   y = f.temp_8_5um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.25 & nbins_y = 60.0
   nchan_used = 2 & channel_wvl_used = [8500,11000]
   tsfc_min = 230.0
   solzen_min = 85.0
end
 
"etropo11_btd8511_day": begin
   classifier_name = 'etropo11_btd8511_day'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'btd8511'
   y = f.temp_8_5um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.25 & nbins_y = 60.0
   nchan_used = 2 & channel_wvl_used = [8500,11000]
   tsfc_min = 230.0
   solzen_max = 95.0
end

"btd1173_btd11133": begin
   classifier_name = 'btd1173_btd11133'
   xtitle = 'btd1173'
   x = f.temp_11_0um_nom - f.temp_7_3um_nom
   x_min = -10.0 & x_bin = 2.0 & nbins_x = 45
   ytitle = 'btd11133'
   y = f.temp_11_0um_nom - f.temp_13_3um_nom
   y_min = -10.0 & y_bin = 2.0 & nbins_y = 40
   nchan_used = 3 & channel_wvl_used = [7300,11000,13300]
end

;--- Begin Sam
"bt67_btd12133": begin
   classifier_name = 'bt67_btd12133'
   xtitle = 'bt67'
   x = f.temp_6_7um_nom
   x_min = 190.0 & x_bin = 5.0 & nbins_x = 20
   ytitle = 'btd12133'
   y = f.temp_12_0um_nom - f.temp_13_3um_nom
   y_min = -10 & y_bin = 2.0 & nbins_y = 25.0
   nchan_used = 3 & channel_wvl_used = [6700,12000,13300]
end

"bt67_btd85133": begin
   classifier_name = 'bt67_btd85133'
   xtitle = 'bt67'
   x = f.temp_6_7um_nom
   x_min = 190.0 & x_bin = 5.0 & nbins_x = 20
   ytitle = 'btd85133'
   y = f.temp_8_5um_nom - f.temp_13_3um_nom
   y_min = -10 & y_bin = 2.0 & nbins_y = 25.0
   nchan_used = 3 & channel_wvl_used = [6700,8500,13300]
end

"bt67_btd1112": begin
   classifier_name = 'bt67_btd1112'
   xtitle = 'bt67'
   x = f.temp_6_7um_nom
   x_min = 190.0 & x_bin = 5.0 & nbins_x = 20
   ytitle = 'btd1112'
   y = f.temp_11_0um_nom - f.temp_12_0um_nom
   y_min = -2 & y_bin = 0.25 & nbins_y = 40.0
   nchan_used = 3 & channel_wvl_used = [6700,11000,12000]
end

"bt67_btd3811_day": begin
   classifier_name = 'bt67_btd3811_day'
   xtitle = 'bt67'
   x = f.temp_6_7um_nom
   x_min = 190.0 & x_bin = 5.0 & nbins_x = 20
   ytitle = 'btd3811'
   y = f.temp_3_75um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.5 & nbins_y = 60.0
   nchan_used = 3 & channel_wvl_used = [3750,6700,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solscatang_min = 95.0
   solglint_mask_min = 0
   solglint_mask_max = 0
   tsfc_min = 240.0
end

"bt67_btd3812_day": begin
   classifier_name = 'bt67_btd3812_day'
   xtitle = 'bt67'
   x = f.temp_6_7um_nom
   x_min = 190.0 & x_bin = 5.0 & nbins_x = 20
   ytitle = 'btd3812'
   y = f.temp_3_75um_nom - f.temp_12_0um_nom
   y_min = -5 & y_bin = 0.5 & nbins_y = 60.0
   nchan_used = 3 & channel_wvl_used = [3750,6700,12000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solscatang_min = 95.0
   solglint_mask_min = 0
   solglint_mask_max = 0
   tsfc_min = 240.0
end

"bt73_btd3885_day": begin
   classifier_name = 'bt73_btd3885_day'
   xtitle = 'bt73'
   x = f.temp_7_3um_nom
   x_min = 190.0 & x_bin = 5.0 & nbins_x = 20
   ytitle = 'btd3885'
   y = f.temp_3_75um_nom - f.temp_8_5um_nom
   y_min = -5 & y_bin = 1.0 & nbins_y = 40.0
   nchan_used = 3 & channel_wvl_used = [3750,7300,8500]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solscatang_min = 95.0
   solglint_mask_min = 0
   solglint_mask_max = 0
   tsfc_min = 240.0
end

"bt73_btd8512": begin
   classifier_name = 'bt73_btd8512'
   xtitle = 'bt73'
   x = f.temp_7_3um_nom
   x_min = 190.0 & x_bin = 5.0 & nbins_x = 20
   ytitle = 'btd8512'
   y = f.temp_8_5um_nom - f.temp_12_0um_nom
   y_min = -10 & y_bin = 0.5 & nbins_y = 40.0
   nchan_used = 3 & channel_wvl_used = [3750,8500,12000]
end

;--- End Sam

"bt11_btd1173": begin
   classifier_name = 'bt11_btd1173'
   xtitle = 'bt11'
   x = f.temp_11_0um_nom
   x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'btd1173'
   y = f.temp_11_0um_nom - f.temp_7_3um_nom
   y_min = -10 & y_bin = 2.0 & nbins_y = 30.0
   nchan_used = 2 & channel_wvl_used = [7300,11000]
end

"bt11_btd1173_day": begin
   classifier_name = 'bt11_btd1173_day'
   xtitle = 'bt11'
   x = f.temp_11_0um_nom
   x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'btd1173'
   y = f.temp_11_0um_nom - f.temp_7_3um_nom
   y_min = -10 & y_bin = 2.0 & nbins_y = 30.0
   nchan_used = 2 & channel_wvl_used = [7300,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
end
 
"bt11_btd1173_night": begin
   classifier_name = 'bt11_btd1173_night'
   xtitle = 'bt11'
   x = f.temp_11_0um_nom
   x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'btd1173'
   y = f.temp_11_0um_nom - f.temp_7_3um_nom
   y_min = -10 & y_bin = 2.0 & nbins_y = 30.0
   nchan_used = 2 & channel_wvl_used = [7300,11000]
   solzen_min = 90.0
   solzen_max = 180
end
 
"bt11_btd11133": begin
   classifier_name = 'bt11_btd11133'
   xtitle = 'bt11'
   x = f.temp_11_0um_nom
   x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'btd11133'
   y = f.temp_11_0um_nom - f.temp_13_3um_nom
   y_min = -10 & y_bin = 2.5 & nbins_y = 30.0
   nchan_used = 2 & channel_wvl_used = [13300,11000]
end
  
"bt11_btd11133_day": begin
   classifier_name = 'bt11_btd11133_day'
   xtitle = 'bt11'
   x = f.temp_11_0um_nom
   x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'btd11133'
   y = f.temp_11_0um_nom - f.temp_13_3um_nom
   y_min = -10 & y_bin = 2.5 & nbins_y = 30.0
   nchan_used = 2 & channel_wvl_used = [13300,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
end

"bt11_btd11133_night": begin
   classifier_name = 'bt11_btd11133_night'
   xtitle = 'bt11'
   x = f.temp_11_0um_nom
   x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'btd11133'
   y = f.temp_11_0um_nom - f.temp_13_3um_nom
   y_min = -10 & y_bin = 2.5 & nbins_y = 30.0
   nchan_used = 2 & channel_wvl_used = [13300,11000]
   solzen_min = 90.0
   solzen_max = 180.0
end

"bt11_btd13373": begin
   classifier_name = 'bt11_btd13373'
   xtitle = 'bt11'
   x = f.temp_11_0um_nom
   x_min = 180.0 & x_bin = 2.0 & nbins_x = 70
   ytitle = 'btd13373'
   y = f.temp_13_3um_nom - f.temp_7_3um_nom
   y_min = -10.0 & y_bin = 1.0 & nbins_y = 30
   nchan_used = 3 & channel_wvl_used = [7300,11000,13300]
end

"btd1112_btd8511": begin
   classifier_name = 'btd1112_btd8511'
   xtitle = 'btd1112'
   x = f.temp_11_0um_nom - f.temp_12_0um_nom
   x_min = -2 & x_bin = 0.25 & nbins_x = 40.0
   ytitle = 'btd8511'
   y = f.temp_8_5um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.25 & nbins_y = 60.0
   nchan_used = 3 & channel_wvl_used = [8500,11000,12000]
   tsfc_min = 240.0
end

"topa_btd8511": begin
   classifier_name = 'topa_btd8511'
   xtitle = 'topa'
   x = f.cld_temp_opaque
   x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'btd8511'
   y = f.temp_8_5um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.25 & nbins_y = 60.0
   nchan_used = 2 & channel_wvl_used = [8500,11000]
   tsfc_min = 230.0
end

"etropo11_btd3811_all": begin
   classifier_name = 'etropo11_btd3811_all'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'btd3811'
   y = f.temp_3_75um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.5 & nbins_y = 60.0
   nchan_used = 2 & channel_wvl_used = [3750,11000]
   tsfc_min = 240.0
   tsfc_min = 230.0
end

"etropo11_btd3811_day": begin
   classifier_name = 'etropo11_btd3811_day'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'btd3811'
   y = f.temp_3_75um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.5 & nbins_y = 60.0
   nchan_used = 2 & channel_wvl_used = [3750,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solscatang_min = 95.0
   solglint_mask_min = 0
   solglint_mask_max = 0
   tsfc_min = 240.0
end

"etropo11_btd1167": begin
   classifier_name = 'etropo11_btd1167'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'btd1167'
   y = f.temp_11_0um_nom - f.temp_6_7um_nom
   y_min = -5 & y_bin = 2.0 & nbins_y = 40.0
   nchan_used = 2 & channel_wvl_used = [6700,11000]
   tsfc_min = 230.0
end

"bt11_btd1167": begin
   classifier_name = 'bt11_btd1167'
   xtitle = 'bt11'
   x = f.temp_11_0um_nom
   x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'btd1167'
   y = f.temp_11_0um_nom - f.temp_6_7um_nom
   y_min = -5 & y_bin = 2.0 & nbins_y = 40.0
   nchan_used = 2 & channel_wvl_used = [6700,11000]
end

"etropo11_btd3811_night": begin
   classifier_name = 'etropo11_btd3811_night'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'btd3811'
   y = f.temp_3_75um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.5 & nbins_y = 60.0
   nchan_used = 2 & channel_wvl_used = [3750,11000]
   solzen_min = 90.0
   solzen_max = 180.0
   solglintzen_min = 0.0
   solglintzen_max = 180.0
   tsfc_min = 240.0
   tsfc_min = 230.0
end

"bt11_btd3811_all": begin
   classifier_name = 'bt11_btd3811_all'
   xtitle = 'bt11'
   x = f.temp_11_0um_nom
   x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'btd3811'
   y = f.temp_3_75um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.5 & nbins_y = 60.0
   nchan_used = 2 & channel_wvl_used = [3750,11000]
   tsfc_min = 240.0
end

"bt11_btd3811_day": begin
   classifier_name = 'bt11_btd3811_day'
   xtitle = 'bt11'
   x = f.temp_11_0um_nom
   x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'btd3811'
   y = f.temp_3_75um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.5 & nbins_y = 60.0
   nchan_used = 2 & channel_wvl_used = [3750,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solscatang_min = 95.0
   solglint_mask_min = 0
   solglint_mask_max = 0
   tsfc_min = 240.0
end

"bt11_btd3811_night": begin
   classifier_name = 'bt11_btd3811_night'
   xtitle = 'bt11'
   x = f.temp_11_0um_nom
   x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'btd3811'
   y = f.temp_3_75um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.5 & nbins_y = 60.0
   nchan_used = 2 & channel_wvl_used = [3750,11000]
   solzen_min = 90.0
   solzen_max = 180.0
   tsfc_min = 240.0
end

"bt11_btd1112": begin
   classifier_name = 'bt11_btd1112'
   xtitle = 'bt11' & x = f.temp_11_0um_nom & x_min = 180 & x_bin = 2.0 & nbins_x = 70
   ytitle = 'btd1112' & y_min = -2 & y_bin = 0.25 & nbins_y = 40.0
   y = f.temp_11_0um_nom - f.temp_12_0um_nom
   nchan_used = 2 & channel_wvl_used = [11000,12000]
   tsfc_min = 240.0
end

"etropo11_btd1112": begin
   classifier_name = 'etropo11_btd1112'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.4 & x_bin = 0.025 & nbins_x = 60
   ytitle = 'btd1112'
   y = f.temp_11_0um_nom - f.temp_12_0um_nom
   y_min = -2 & y_bin = 0.25 & nbins_y = 40.0
   nchan_used = 2 & channel_wvl_used = [11000,12000]
   tsfc_min = 240.0
end

"topa_btd1112": begin
   classifier_name = 'topa_btd1112'
   xtitle = 'topa'
   x = f.cld_temp_opaque
   x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'btd1112'
   y = f.temp_11_0um_nom - f.temp_12_0um_nom
   y_min = -2 & y_bin = 0.25 & nbins_y = 40.0
   nchan_used = 2 & channel_wvl_used = [11000,12000]
   tsfc_min = 240.0
end

"btdclr11_logbt11std": begin
   classifier_name ='btdclr11_logbt11std'
   xtitle = 'btdclr11'
   x = f.temp_11_0um_nom_clear_sky - f.temp_11_0um_nom
   x_min = -15 & x_bin = 0.5 & nbins_x = 100.0
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [11000]
   coast_mask_max = 0
;  zsfc_max = 2000.0
   zsfc_std_max = zsfc_std_max_iruni
end

"btdclr11_logbt11std_night": begin
   classifier_name ='btdclr11_logbt11std_night'
   xtitle = 'btdclr11'
   x = f.temp_11_0um_nom_clear_sky - f.temp_11_0um_nom
   x_min = -15 & x_bin = 0.5 & nbins_x = 100.0
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [11000]
   coast_mask_max = 0
;  zsfc_max = 2000.0
   zsfc_std_max = zsfc_std_max_iruni
   solzen_min = 90.0
end
"btdclr11_logbt11std_day": begin
   classifier_name ='btdclr11_logbt11std_day'
   xtitle = 'btdclr11'
   x = f.temp_11_0um_nom_clear_sky - f.temp_11_0um_nom
   x_min = -15 & x_bin = 0.5 & nbins_x = 100.0
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [11000]
   coast_mask_max = 0
;  zsfc_max = 2000.0
   zsfc_std_max = zsfc_std_max_iruni
   solzen_max = 90.0
end

"bt11_logbt11std": begin
   classifier_name ='bt11_logbt11std'
   xtitle = 'bt11'
   x = f.temp_11_0um_nom
   x_min = 180.0 & x_bin = 2.0 & nbins_x = 80
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [11000]
   coast_mask_max = 0
;  zsfc_max = 2000.0
   zsfc_std_max = zsfc_std_max_iruni
end

"bt10_logbt10std": begin
   classifier_name ='bt10_logbt10std'
   xtitle = 'bt10'
   x = f.temp_10_4um_nom
   x_min = 180.0 & x_bin = 2.0 & nbins_x = 80
   ytitle = 'logbt10std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_10_4UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [10400]
   coast_mask_max = 0
;  zsfc_max = 2000.0
   zsfc_std_max = zsfc_std_max_iruni
   coast_mask_max = 0
end

"etropo11_logbt11std": begin
   classifier_name ='etropo11_logbt11std'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.6 & x_bin = 0.025 & nbins_x = 80
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [11000]
   coast_mask_max = 0
   tsfc_min = 230.0
   zsfc_std_max = zsfc_std_max_iruni
end

"etropo11_logbt11std_day": begin
   classifier_name ='etropo11_logbt11std_day'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.6 & x_bin = 0.025 & nbins_x = 80
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [11000]
   solzen_min=0.0
   ;solzen_max=95.0
   coast_mask_max = 0
   tsfc_min = 230.0
   zsfc_std_max = zsfc_std_max_iruni
end

"etropo11_logbt11std_night": begin
   classifier_name ='etropo11_logbt11std_night'
   xtitle = 'etropo11'
   x = f.emiss_tropo_11_0um_nom
   x_min = -0.6 & x_bin = 0.025 & nbins_x = 80
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [11000]
   solzen_min=85.0
   solzen_max=180.0
   coast_mask_max = 0
   tsfc_min = 230.0
   zsfc_std_max = zsfc_std_max_iruni
end

"etropo10_logbt10std": begin
   classifier_name ='etropo10_logbt10std'
   xtitle = 'etropo10'
   x = f.emiss_tropo_10_4um_nom
   x_min = -0.6 & x_bin = 0.025 & nbins_x = 80
   ytitle = 'logbt10std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_10_4UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [10400]
   coast_mask_max = 0
   tsfc_min = 230.0
   zsfc_std_max = zsfc_std_max_iruni
end

"topa_logbt11std": begin
   classifier_name ='topa_logbt11std'
   xtitle = 'topa'
   x = f.cld_temp_opaque
   x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'logbt11std'
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   y_min = -1.5 & y_bin = 0.25 & nbins_y = 15.0
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [11000]
   coast_mask_max = 0
   zsfc_std_max = zsfc_std_max_iruni
end

"refl160_refl138": begin
   classifier_name = 'refl160_refl138'
   xtitle = 'refl160'
   x = f.refl_1_60um_nom
   x_min = -2.0 & x_bin = 1.0 & nbins_x = 65
   ytitle = 'refl138'
   y = f.refl_1_38um_nom
   y_min = -2.0 & y_bin = 0.5 & nbins_y = 80.0
   nchan_used = 2 & channel_wvl_used = [1380,1600]
   solglint_mask_min = 0
   solglint_mask_max = 0
   slant_tpw_min = 0.5
end

"tsfc_btd1112": begin
   classifier_name = 'tsfc_btd1112'
   xtitle = 'tsfc'
   x = f.surface_temperature_nwp
   x_min = 200.0 & x_bin = 5.0 & nbins_x = 25
   ytitle = 'btd1112'
   y = f.temp_11_0um_nom - f.temp_12_0um_nom
   y_min = -2.0 & y_bin = 0.2 & nbins_y = 50
   nchan_used = 2 & channel_wvl_used = [11000,12000]
   tsfc_min = 240.0
end

"dtsfcbt11_btd1112": begin
   classifier_name = 'dtsfcbt11_btd1112'
   xtitle = 'dtsfcbt11'
   x = f.surface_temperature_nwp - f.temp_11_0um_nom
   x_min = -20.0 & x_bin = 5.0 & nbins_x = 30
   ytitle = 'btd1112'
   y = f.temp_11_0um_nom - f.temp_12_0um_nom
   y_min = -2.0 & y_bin = 0.2 & nbins_y = 50
   nchan_used = 2 & channel_wvl_used = [11000,12000]
   tsfc_min = 240.0
end

"zopa_btd1112": begin
   classifier_name = 'zopa_btd1112'
   xtitle = 'zopa'
   x = f.cld_height_opaque
   x_min = 0.0 & x_bin = 1000.0 & nbins_x = 20
   ytitle = 'btd1112'
   y = f.temp_11_0um_nom - f.temp_12_0um_nom
   y_min = -2.0 & y_bin = 0.2 & nbins_y = 50
   nchan_used = 2 & channel_wvl_used = [11000,12000]
   tsfc_min = 240.0
end

"zopa_fmft": begin
   classifier_name = 'zopa_fmft'
   xtitle = 'zopa'
   x = f.cld_height_opaque
   x_min = 0.0 & x_bin = 1000.0 & nbins_x = 20
   ytitle = 'fmft'
   y = fmft(f.temp_11_0um_nom, f.temp_12_0um_nom, f.temp_11_0um_nom_clear_sky, f.temp_12_0um_nom_clear_sky)
   y_min = -2.0 & y_bin = 0.2 & nbins_y = 50
   nchan_used = 2 & channel_wvl_used = [11000,12000]
   tsfc_min = 240.0
end

"zopa_btd8511": begin
   classifier_name = 'zopa_btd8511'
   xtitle = 'zopa'
   x = f.cld_height_opaque
   x_min = 0.0 & x_bin = 1000.0 & nbins_x = 20
   ytitle = 'btd8511'
   y = f.temp_8_5um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.25 & nbins_y = 60.0
   nchan_used = 2 & channel_wvl_used = [8500,11000]
   tsfc_min = 230.0
end

"logzopa_btd8511": begin
   classifier_name = 'logzopa_btd8511'
   xtitle = 'logzopa'
   xtitle = 'logzopa' & x = f.cld_height_opaque & x_min = 1.5 & x_bin = 0.1 & nbins_x = 30.0
   idx = where(x gt 0.0,cc) & if (cc gt 0) then x[idx] = alog10(x[idx])
   idx = where(x eq 0.0,cc) & if (cc gt 0) then x[idx] = x_min
   ytitle = 'btd8511'
   y = f.temp_8_5um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.25 & nbins_y = 60.0
   nchan_used = 2 & channel_wvl_used = [8500,11000]
   tsfc_min = 230.0
end

"bt11_btd8511": begin
   classifier_name = 'bt11_btd8511'
   xtitle = 'bt11_btd8511'
   xtitle = 'bt11'
   x = f.temp_11_0um_nom
   x_min = 180.0 & x_bin = 5.0 & nbins_x = 32
   ytitle = 'btd8511'
   y = f.temp_8_5um_nom - f.temp_11_0um_nom
   y_min = -5 & y_bin = 0.25 & nbins_y = 60.0
   nchan_used = 2 & channel_wvl_used = [8500,11000]
   tsfc_min = 230.0
end

"btdclr11_logbt11std": begin
   classifier_name = 'btdclr11_logbt11std'
   xtitle = 'btdclr11'
   x = f.temp_11_0um_nom_clear_sky - f.temp_11_0um_nom
   x_min = -15 & x_bin = 0.5 & nbins_x = 100.0
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [11000]
;  zsfc_max = 2000.0
   zsfc_std_max = zsfc_std_max_iruni
   coast_mask_max = 0
end

"btdclr11_logbt11std_night": begin
   classifier_name = 'btdclr11_logbt11std_night'
   xtitle = 'btdclr11'
   x = f.temp_11_0um_nom_clear_sky - f.temp_11_0um_nom
   x_min = -15 & x_bin = 0.5 & nbins_x = 100.0
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [11000]
;  zsfc_max = 2000.0
   zsfc_std_max = zsfc_std_max_iruni
   coast_mask_max = 0
   solzen_min = 90.0
end

"logzopa_logbt11std": begin
   classifier_name = 'logzopa_logbt11std'
   xtitle = 'logzopa' & x = f.cld_height_opaque & x_min = 1.5 & x_bin = 0.1 & nbins_x = 30.0
   idx = where(x gt 0.0,cc) & if (cc gt 0) then x[idx] = alog10(x[idx])
   idx = where(x eq 0.0,cc) & if (cc gt 0) then x[idx] = x_min
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [11000]
   zsfc_std_max = zsfc_std_max_iruni
   coast_mask_min = 0
   coast_mask_max = 0
end

"logzopa_logdtsfcopa": begin
   classifier_name = 'logzopa_logdtsfcopa'
   xtitle = 'logzopa' & x = f.cld_height_opaque & x_min = 1.5 & x_bin = 0.1 & nbins_x = 30.0
   idx = where(x gt 0.0,cc) & if (cc gt 0) then x[idx] = alog10(x[idx])
   idx = where(x eq 0.0,cc) & if (cc gt 0) then x[idx] = x_min
   ytitle = 'logdtsfcopa'
   y_min = 0.1 & y_bin = 0.1 & nbins_y = 20.0
   y = f.surface_temperature_nwp - f.cld_temp_opaque
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [11000]
   zsfc_std_max = zsfc_std_max_iruni
   coast_mask_min = 0
   coast_mask_max = 0
end

"logzopa_dtsfcopa": begin
   classifier_name = 'logzopa_dtsfcopa'
   xtitle = 'logzopa' & x = f.cld_height_opaque & x_min = 1.5 & x_bin = 0.1 & nbins_x = 30.0
   idx = where(x gt 0.0,cc) & if (cc gt 0) then x[idx] = alog10(x[idx])
   idx = where(x eq 0.0,cc) & if (cc gt 0) then x[idx] = x_min
   ytitle = 'dtsfcopa'
   y_min = -10.0 & y_bin = 2.0 & nbins_y = 30.0
   y = f.surface_temperature_nwp - f.cld_temp_opaque
   nchan_used = 1 & channel_wvl_used = [11000]
   zsfc_std_max = zsfc_std_max_iruni
   coast_mask_min = 0
   coast_mask_max = 0
end

"zopa_logbt10std": begin
   classifier_name = 'zopa_logbt10std'
   xtitle = 'zopa'
   x = f.cld_height_opaque
   x_min = 0.0 & x_bin = 1000.0 & nbins_x = 20
   ytitle = 'logbt10std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_10_4UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 1 & channel_wvl_used = [10400]
   coast_mask_max = 0
   zsfc_std_max = zsfc_std_max_iruni
   coast_mask_max = 0
end

"ndvi_ndsi": begin
   classifier_name = 'ndvi_ndsi'
   xtitle = 'ndvi'
   x = ndvi(f.refl_0_65um_nom, f.refl_0_86um_nom, missing) & x_min = -1.0 & x_bin = 0.05 & nbins_x = 50

   ytitle = 'ndsi'
   y = ndsi(f.refl_0_65um_nom, f.refl_1_60um_nom, missing) & y_min = -1.0 & y_bin = 0.05 & nbins_y = 50

   nchan_used = 3 & channel_wvl_used = [650,860,1600]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solscatang_min = 95.0
   solglint_mask_min = 0
   solglint_mask_max = 0
   coast_mask_max = 0
end

"drefl065clr_ndsi": begin
   classifier_name = 'drefl065clr_ndsi'
   xtitle = 'drefl065clr'
   x = f.refl_0_65um_nom - f.refl_0_65um_nom_clear_sky & x_min = -40.0 & x_bin = 1.0 & nbins_x = 120

   ytitle = 'ndsi'
   y = ndsi(f.refl_0_65um_nom, f.refl_1_60um_nom, missing) & y_min = -0.5 & y_bin = 0.025 & nbins_y = 60

   nchan_used = 2 & channel_wvl_used = [650,1600]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solscatang_min = 95.0
   solglint_mask_min = 0
   solglint_mask_max = 0
   coast_mask_max = 0
   snow_class_max = 1
end

"refl065_drefl065clr": begin
    classifier_name = 'refl065_drefl065clr'
    xtitle = 'refl065'
    x = f.refl_0_65um_nom & x_min = -10.0 & x_bin = 2.0 & nbins_x = 50
    ytitle = 'drefl065clr'
    y = f.refl_0_65um_nom - f.refl_0_65um_nom_clear_sky  & y_min = -20.0 & y_bin = 2.0 & nbins_y = 40.0
    nchan_used = 1 & channel_wvl_used = [650]
    solzen_min = 0.0
    solzen_max = solzen_day_max
    solglint_mask_min = 0
    solglint_mask_max = 0
    snow_class_max = 1
end

"refl047_refl160": begin
    classifier_name = 'refl047_refl160'
    xtitle = 'refl047'
    x = f.refl_0_47um_nom & x_min = -10.0 & x_bin = 2.0 & nbins_x = 50
    ytitle = 'refl160'
    y = f.refl_1_60um_nom & y_min = -10.0 & y_bin = 2.0 & nbins_y = 50
    nchan_used = 2 & channel_wvl_used = [650,1600]
    solzen_min = 0.0
    solzen_max = solzen_day_max
    solglint_mask_min = 0
    solglint_mask_max = 0
    snow_class_max = 1
end

"refl065_refl065cv": begin
   classifier_name = 'refl065_refl065cv'
   xtitle = 'refl065'
   x = f.refl_0_65um_nom & x_min = -10.0 & x_bin = 2.0 & nbins_x = 40
   ytitle = 'refl065cv'
   y = f.refl_0_65um_nom_stddev_sub / f.refl_0_65um_nom
   idx = where(f.refl_0_65um_nom_stddev_sub eq missing or f.refl_0_65um_nom eq missing, cc)
   if (cc gt 0) then y[idx] = missing
   y_min = 0.0 & y_bin = 0.02 & nbins_y = 50.0
   nchan_used = 1 & channel_wvl_used = [650]
   coast_mask_max = 0
   solzen_min = 0.0
   solzen_max = solzen_day_max
   snow_class_max = 1
;  zsfc_max = 2000.0
  end

"refl065_logrefl065std": begin
    classifier_name ='refl065_logrefl065std'
    xtitle = 'refl065'
    x = f.refl_0_65um_nom & x_min = -10.0 & x_bin = 2.0 & nbins_x = 50
    ytitle = 'logrefl065std'
    y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
    y = f.REFL_0_65UM_NOM_STDDEV_3X3
    idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
    idx = where(y eq 0.0,cc) 
    if (cc gt 0) then y[idx] = y_min
    nchan_used = 1 & channel_wvl_used = [650]
    coast_mask_min = 0
    coast_mask_max = 0
    solzen_min = 0.0
    solzen_max = solzen_day_max
    snow_class_max = 1
;   zsfc_max = 2000.0
end

"refl160_logrefl065std": begin
    classifier_name ='refl160_logrefl065std'
    xtitle = 'refl160'
    x = f.refl_1_60um_nom & x_min = -10.0 & x_bin = 2.0 & nbins_x = 50
    ytitle = 'logrefl065std'
    y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
    y = f.REFL_0_65UM_NOM_STDDEV_3X3
    idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
    idx = where(y eq 0.0,cc) 
    if (cc gt 0) then y[idx] = y_min
    nchan_used = 2 & channel_wvl_used = [650,1600]
    coast_mask_min = 0
    coast_mask_max = 0
    solzen_min = 0.0
    solzen_max = solzen_day_max
;   zsfc_max = 2000.0
end

"drefl065clr_refl065std": begin
    classifier_name = 'drefl065clr_refl065std'
    xtitle = 'drefl065clr' & x = f.refl_0_65um_nom - f.refl_0_65um_nom_clear_sky & x_min = -40.0 & x_bin = 1.0 & nbins_x = 120
    ytitle = 'refl065std' & y_min = 0.0 & y_bin = 0.10 & nbins_y = 100.0 & y = f.REFL_0_65UM_NOM_STDDEV_3X3
    nchan_used = 2 & channel_wvl_used = [650,11000]
    coast_mask_min = 0
    coast_mask_max = 0
    solzen_min = 0.0
    solzen_max = solzen_day_max
    snow_class_max = 1
;   zsfc_max = 2000.0
end

"drefl065clr_logrefl065std": begin
    classifier_name = 'drefl065clr_logrefl065std'
    xtitle = 'drefl065clr'
    x = f.refl_0_65um_nom - f.refl_0_65um_nom_clear_sky & x_min = -40.0 & x_bin = 1.0 & nbins_x = 120
    ytitle = 'logrefl065std'
    y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
    y = f.REFL_0_65UM_NOM_STDDEV_3X3
    idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
    idx = where(y eq 0.0,cc) 
    if (cc gt 0) then y[idx] = y_min
    nchan_used = 2 & channel_wvl_used = [650,11000]
    coast_mask_min = 0
    coast_mask_max = 0
    solzen_min = 0.0
    solzen_max = solzen_day_max
    snow_class_max = 1
;   zsfc_max = 2000.0
end

"logbt11std_logrefl065std": begin
    classifier_name ='logbt11std_logrefl065std'
    xtitle = 'logbt11std'
    ;x_min = -1.5 & x_bin = 0.10 & nbins_x = 36.0
    x_min = -2.0 & x_bin = 0.10 & nbins_x = 40.0
    x = f.TEMP_11_0UM_NOM_STDDEV_3X3
    idx = where(x gt 0.0,cc) & if (cc gt 0) then x[idx] = alog10(x[idx])
    idx = where(x eq 0.0,cc) & if (cc gt 0) then x[idx] = x_min
    ytitle = 'logrefl065std'
    ;y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
    y_min = -2.0 & y_bin = 0.10 & nbins_y = 40.0
    y = f.REFL_0_65UM_NOM_STDDEV_3X3
    idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
    idx = where(y eq 0.0,cc) 
    if (cc gt 0) then y[idx] = y_min
    nchan_used = 2 & channel_wvl_used = [650,11000]
    coast_mask_min = 0
    coast_mask_max = 0
    solzen_min = 0.0
    solzen_max = solzen_day_max
  ; zsfc_max = 2000.0
    zsfc_std_max = zsfc_std_max_iruni
end

"bt11std_refl065std": begin
    classifier_name ='bt11std_refl065std'
    xtitle = 'bt11std' & x = f.temp_11_0um_nom_stddev_3x3 & x_min = 0.0 & x_bin = 0.1 & nbins_x = 100.0
    ytitle = 'refl065std' & y = f.refl_0_65um_nom_stddev_3x3 & y_min = 0.0 & y_bin = 0.1 & nbins_y = 100.0
    nchan_used = 2 & channel_wvl_used = [650,11000]
    coast_mask_min = 0
    coast_mask_max = 0
    solzen_min = 0.0
    solzen_max = solzen_day_max
;   zsfc_max = 2000.0
    zsfc_std_max = zsfc_std_max_iruni
end

"btd3811_logbt11std_all": begin
   classifier_name = 'btd3811_logbt11std_all'
   xtitle = 'btd3811'
   x = f.temp_3_75um_nom - f.temp_11_0um_nom
   x_min = -5 & x_bin = 0.5 & nbins_x = 60.0
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 2 & channel_wvl_used = [3750,11000]
;  zsfc_max = 2000.0
   zsfc_std_max = zsfc_std_max_iruni
   tsfc_min = 240.0
end

"btd3811_logbt11std_day": begin
   classifier_name = 'btd3811_logbt11std_day'
   xtitle = 'btd3811'
   x = f.temp_3_75um_nom - f.temp_11_0um_nom
   x_min = -5 & x_bin = 0.5 & nbins_x = 60.0
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 2 & channel_wvl_used = [3750,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   solscatang_min = 95.0
   solglint_mask_min = 0
   solglint_mask_max = 0
;  zsfc_max = 2000.0
   zsfc_std_max = zsfc_std_max_iruni
   tsfc_min = 240.0
end

"btd3811_logbt11std_night": begin
   classifier_name = 'btd3811_logbt11std_night'
   xtitle = 'btd3811'
   x = f.temp_3_75um_nom - f.temp_11_0um_nom
   x_min = -5 & x_bin = 0.5 & nbins_x = 60.0
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 2 & channel_wvl_used = [3750,11000]
   solzen_min = 90.0
   solzen_max = 180.0
;  zsfc_max = 2000.0
   zsfc_std_max = zsfc_std_max_iruni
   tsfc_min = 240.0
end

"dbtd3811clr_logbt11std_day": begin
   classifier_name = 'dbtd3811clr_day'
   xtitle = 'dbtd3811clr'
   x = (f.temp_3_75um_nom - f.temp_11_0um_nom) - $
       (f.temp_3_75um_nom_clear_sky - f.temp_11_0um_nom_clear_sky)
   x_min = -10 & x_bin = 0.5 & nbins_x = 80.0
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 2 & channel_wvl_used = [3750,11000]
   solzen_min = 0.0
   solzen_max = solzen_day_max
   tsfc_min = 240.0
   zsfc_std_max = zsfc_std_max_iruni
  end

"dbtd3811clr_logbt11std_night": begin
   classifier_name = 'dbtd3811clr_night'
   xtitle = 'dbtd3811clr'
   x = (f.temp_3_75um_nom - f.temp_11_0um_nom) - $
       (f.temp_3_75um_nom_clear_sky - f.temp_11_0um_nom_clear_sky)
   x_min = -10 & x_bin = 0.5 & nbins_x = 80.0
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 2 & channel_wvl_used = [3750,11000]
   solzen_min = 90.0
   solzen_max = 180.0
   tsfc_min = 240.0
   zsfc_std_max = zsfc_std_max_iruni
  end

"bt6711covar_logbt11std": begin
   classifier_name = 'bt6711covar'
   xtitle = 'bt6711covar'
   x = f.temp_11um_vs_67um_covar_5x5
   x_min = -4 & x_bin = 0.15 & nbins_x = 80
   ytitle = 'logbt11std'
   y_min = -1.5 & y_bin = 0.10 & nbins_y = 36.0
   y = f.TEMP_11_0UM_NOM_STDDEV_3X3
   idx = where(y gt 0.0,cc) & if (cc gt 0) then y[idx] = alog10(y[idx])
   idx = where(y eq 0.0,cc) & if (cc gt 0) then y[idx] = y_min
   nchan_used = 2 & channel_wvl_used = [6700,11000]
;  zsfc_max = 2000.0
   zsfc_std_max = zsfc_std_max_iruni
end

"bt6711covar_bt11": begin

   classifier_name = 'bt6711covar'
   xtitle = 'bt6711covar'
   x = f.temp_11um_vs_67um_covar_5x5
   x_min = -4 & x_bin = 0.15 & nbins_x = 80
   ytitle = 'bt11'
   y = f.temp_11_0um_nom
   y_min = 180.0 & y_bin = 5.0 & nbins_y = 32
   nchan_used = 2 & channel_wvl_used = [6700,11000]
end

else: begin
      print, "error: unknown test"
      stop
      end

endcase

;------------------------------------------------
; fake other dims 
;------------------------------------------------
ztitle = "blank"
nbins_z = 1
z_min = 0.0
z_bin = 0.0

;------------------------------------------------
; 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
;----------------------------------------------
sfc_type_mask = f.bayes_mask_sfc_type

if (keyword_set(isccp)) then begin
   make_ecm_surface_type, f, sfc_type_mask_2, /isccp
endif


;-----------------------------------------------
;  apply filters
;----------------------------------------------
cod_lidar = f.CLOSEST_CALIPSO_COD
slant_tpw = f.total_precipitable_water_nwp / cos(f.sensor_zenith_angle * !DTOR)

if (keyword_set(new)) then begin
      frac_lidar = f.CLOSEST_CALIPSO_CLOUD_FRACTION_5KM
      tc_lidar = f.CLOSEST_CALIPSO_TOP_MID_TEMPERATURE
      phase_lidar = f.CLOSEST_CALIPSO_CLOUD_PHASE_TOP_LAYER
      ;frac_lidar = f.CLOSEST_CALIPSO_CLOUD_FRACTION
      ;tc_lidar = f.CLOSEST_CALIPSO_MID_TEMPERATURE
      ;phase_lidar = f.CLOSEST_CALIPSO_PHASE_TOP_LAYER
      codi_lidar = f.CLOSEST_CALIPSO_COD_ICE
endif 
if (keyword_set(old)) then begin
      frac_lidar =f.CLOSEST_CALIPSO_CLOUD_FRACTION_x5
      codi_lidar = f.CLOSEST_CALIPSO_COD
      codi_lidar[*] = missing
endif 

;--- convert tc_lidar from C to K
idx = where(tc_lidar ne  missing, cc)
if (cc gt 0) then tc_lidar[idx] = tc_lidar[idx] + 273.15

make_true_phase, phase_lidar, frac_lidar, tc_lidar, cod_lidar, codi_lidar,  $
                    cod_clear_water_thresh, cod_clear_ice_thresh, $
                    f.bayes_mask_sfc_type, f.cld_temp_opaque, f.cld_height_opaque, $
                    f.solar_zenith_angle, f.glint_mask, f.coast_mask, $
                    f.surface_elevation, f.surface_temperature_nwp, $
                    f.refl_0_65um_nom_min_3x3, f.refl_0_65um_nom, $
                    f.temp_11_0um_nom_max_3x3, f.temp_11_0um_nom, $
                    f.cld_opd_mask_0_65um_nom, cod_cloud_thresh, $
                    f.emiss_tropo_11_0um_nom, etrop_cloud_thresh, mod_flag, phase_true, missing

;-------------------------------------------------------------------------------
; fill in things that might be missing in training data
;-------------------------------------------------------------------------------
idx = where(tag_names(f) eq "CITY_MASK",cc)
if (cc eq 0) then begin
   city_mask_temp = make_array(ntrain,/FLOAT, VALUE=MISSING)
endif else begin
   city_mask_temp = f.city_mask
endelse

idx = where(tag_names(f) eq "MOON_ILLUM_FRAC",cc)
if (cc eq 0) then begin
   moon_illum_frac_temp = make_array(ntrain,/FLOAT, VALUE=MISSING)
endif else begin
   moon_illum_frac_temp = f.moon_illum_frac
endelse

apply_filters, phase_true, $
               f.CLOSEST_CALIPSO_TIME_DIFFERENCE, Time_Diff_Max, $
               f.sensor_zenith_angle, zen_min, zen_max, $
               f.solar_zenith_angle, solzen_min, solzen_max, $
               f.lunar_zenith_angle, lunzen_min, lunzen_max, $
               f.land_class,  $
	       frac_lidar, cfrac_min, cfrac_max,  $
               f.glint_zenith_angle, solglintzen_min, solglintzen_max,  $
               f.lunar_glint_zenith_angle, lunglintzen_min, lunglintzen_max, $
               f.scattering_angle, solscatang_min, solscatang_max,  $
               f.snow_class, snow_class_min, snow_class_max, $
               f.glint_mask, solglint_mask_min, solglint_mask_max, $
               f.glint_mask_lunar, lunglint_mask_min, lunglint_mask_max, $
               city_mask_temp, city_mask_min, city_mask_max, $
               moon_illum_frac_temp, moon_illum_frac_min, moon_illum_frac_max, $
               f.coast_mask, coast_mask_min, coast_mask_max, $
               slant_tpw, slant_tpw_min, slant_tpw_max, $
               f.surface_temperature_nwp, tsfc_min, tsfc_max, $
               f.surface_elevation, zsfc_min, zsfc_max, $
               f.surface_elevation_stddev_3x3, zsfc_std_min, zsfc_std_max, $
               ;f.latitude, lat_min, lat_max, missing, /phase_uni_flag
               f.latitude, lat_min, lat_max, missing

;-------------------------------
; select valid data within these limits
;-------------------------------
idx = where(phase_true ge 0.0 and x ne missing and y ne missing,cc)

print, 'total number after filtering= ', cc

;--- throw away filtered data
x = x[idx]
y = y[idx]
phase_true = phase_true[idx]
sfc_type_mask = sfc_type_mask[idx]

;-------------------------------------------------------------
; make arrays
;-------------------------------------------------------------
ice_table = make_array(nbins_x, nbins_y, nsfc, /FLOAT, VALUE=MISSING)
water_table = make_array(nbins_x, nbins_y, nsfc, /FLOAT, VALUE=MISSING)
clear_table = make_array(nbins_x, nbins_y, nsfc, /FLOAT, VALUE=MISSING)
obs_table = make_array(nbins_x, nbins_y, nsfc, /FLOAT, VALUE=0.0)
cloud_fraction = make_array(nsfc,/FLOAT,VALUE=missing)
clear_fraction = make_array(nsfc,/FLOAT,VALUE=missing)
ice_fraction = make_array(nsfc,/FLOAT,VALUE=missing)
water_fraction = make_array(nsfc,/FLOAT,VALUE=missing)
observation_count = make_array(nsfc,/FLOAT,VALUE=missing)

;--------------------------------------------------------------------
; loop over surface type
;--------------------------------------------------------------------
;for isfc = 0, 0 do begin             ;sfc_loop
for isfc = 0, nsfc - 1 do begin             ;sfc_loop

idx = where(sfc_type_mask eq isfc + 1,cc)
if (cc gt count_min) then begin
   x_sfc = x[idx]
   y_sfc = y[idx]
   phase_true_sfc = phase_true[idx]
   observation_count[isfc] = cc
endif else begin
   goto, cycle_sfc
endelse

;---  make priors
idx = where(phase_true_sfc ge 0, count_all)
idx = where(phase_true_sfc eq 2, count_yes_ice_all)
idx = where(phase_true_sfc eq 1, count_yes_water_all)
idx = where(phase_true_sfc eq 0, count_yes_clear_all)

cloud_fraction[isfc] = float(count_yes_water_all + count_yes_ice_all) / float(count_all)
ice_fraction[isfc] = float(count_yes_ice_all) / float(count_all)
water_fraction[isfc] = float(count_yes_water_all) / float(count_all)
clear_fraction[isfc] = 1.0 - cloud_fraction[isfc]

cond_ratio_ice = make_array(nbins_x, nbins_y, /FLOAT, VALUE=MISSING)
cond_ratio_water = make_array(nbins_x, nbins_y, /FLOAT, VALUE=MISSING)
cond_ratio_clear = make_array(nbins_x, nbins_y, /FLOAT, VALUE=MISSING)
post_prob_ice = make_array(nbins_x, nbins_y, /FLOAT, VALUE=MISSING)
post_prob_water = make_array(nbins_x, nbins_y, /FLOAT, VALUE=MISSING)
post_prob_clear = make_array(nbins_x, nbins_y, /FLOAT, VALUE=MISSING)
obs_prob_lidar = make_array(nbins_x, nbins_y, /FLOAT, VALUE=0.0)

for i = 0, nbins_x  -1 do begin     ;i point_loop

    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     ;j-point loop

       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 lt x_max_bin and $
                       y_sfc ge y_min_bin and y_sfc lt y_max_bin and $
                       phase_true_sfc ge 0, count_all_bin)

       if (count_all_bin le 0) then goto, skip_bin
       if (count_all_bin lt count_min) then goto, skip_bin

       idx_ice = where(phase_true_sfc[idx_all] eq 2, count_yes_ice_bin)
       idx_water = where(phase_true_sfc[idx_all] eq 1, count_yes_water_bin)
       idx_clear = where(phase_true_sfc[idx_all] eq 0, count_yes_clear_bin)

       ;--- clear
       count_yes_bin = count_yes_clear_bin
       count_yes_all = count_yes_clear_all
       count_no_bin = count_all_bin - count_yes_clear_bin
       count_no_all = count_all - count_yes_clear_all
       prior_yes = clear_fraction[isfc] 
       compute_class_cond, nb_rank, laplace_flag, $
                        count_yes_all, count_no_all, count_yes_bin, count_no_bin, $
                        max_r, prior_yes, cond_yes, cond_no, r, posterior_prob
       if (r lt 0.0) then begin
	         print, 'negative clear r'
		 stop
       endif

       cond_ratio_clear[i,j] = r
       post_prob_clear[i,j] = posterior_prob

       ;--- water
       count_yes_bin = count_yes_water_bin
       count_yes_all = count_yes_water_all
       count_no_bin = count_all_bin - count_yes_water_bin
       count_no_all = count_all - count_yes_water_all
       prior_yes = water_fraction[isfc] 
       compute_class_cond, nb_rank, laplace_flag, $
                        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_water[i,j] = r
       post_prob_water[i,j] = posterior_prob

       if (r lt 0.0) then begin
	         print, 'negative water r'
		 stop
       endif

       ;--- ice
       count_yes_bin = count_yes_ice_bin
       count_yes_all = count_yes_ice_all
       count_no_bin = count_all_bin - count_yes_ice_bin
       count_no_all = count_all - count_yes_ice_all
       prior_yes = ice_fraction[isfc] 
       compute_class_cond, nb_rank, laplace_flag, $
                        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_ice[i,j] = r
       post_prob_ice[i,j] = posterior_prob

       if (r lt 0.0) then begin
	         print, 'negative ice r'
		 stop
       endif

       if (count_all gt count_min) then begin
        obs_prob_lidar[i,j] = float(count_all_bin)    
       endif

    skip_bin:

    endfor    ;end of j-point loop
endfor    ;end of i-point loop

;--- normalize obs_prob
idx = where(obs_prob_lidar gt 0,cc)
if (cc gt 0) then obs_prob_lidar[idx] = obs_prob_lidar[idx] / total(obs_prob_lidar)

;goto, skip_post_process
obs_prob_thresh = 3.0 / nbins_x / nbins_y

;-------------------------------------------------------------------------
; fill in holes
;-------------------------------------------------------------------------
;goto, skip_fill
obs_prob_interp_value = 1.0e-06
;idx = where(obs_prob_lidar eq 0,cc) & print, 'number with zero before fill = ', cc
;print, 'sum of obs_prob_lidar before fill = ', total(obs_prob_lidar)

;for i = 0, 100  do begin

for ix = 0, nbins_x - 1 do begin

ix_start = max([0,ix-1])
ix_end = min([nbins_x-1,ix+1])

for iy = 0, nbins_y - 1 do begin

    iy_start = max([0,iy-1])
    iy_end = min([nbins_y-1,iy+1])

    if (obs_prob_lidar[ix,iy] eq 0) then begin

       v = obs_prob_lidar[ix_start:ix_end,iy_start:iy_end]
       idx = where(v gt 0,cc)

       if (cc gt 0) then begin

          obs_prob_lidar[ix,iy] = mean(v[idx]) ; or some special number?
          ;obs_prob_lidar[ix,iy] = obs_prob_interp_value ;mean(v[idx]) ; or some special number?

          v = cond_ratio_clear[ix_start:ix_end,iy_start:iy_end]
          cond_ratio_clear[ix,iy] = mean(v[idx])

          v = cond_ratio_water[ix_start:ix_end,iy_start:iy_end]
          cond_ratio_water[ix,iy] = mean(v[idx])

          v = cond_ratio_ice[ix_start:ix_end,iy_start:iy_end]
          cond_ratio_ice[ix,iy] = mean(v[idx])

          v = post_prob_clear[ix_start:ix_end,iy_start:iy_end]
          post_prob_clear[ix,iy] = mean(v[idx])

          v = post_prob_water[ix_start:ix_end,iy_start:iy_end]
          post_prob_water[ix,iy] = mean(v[idx])

          v = post_prob_ice[ix_start:ix_end,iy_start:iy_end]
          post_prob_ice[ix,iy] = mean(v[idx])

       endif
    endif

endfor
endfor
;endfor
skip_fill:

;idx = where(obs_prob_lidar eq 0,cc) & print, 'number with zero after fill = ', cc
;print, 'sum of obs_prob_lidar after fill = ', total(obs_prob_lidar)

;--- renormalize obs_prob
idx = where(obs_prob_lidar gt 0,cc)
if (cc gt 0) then obs_prob_lidar[idx] = obs_prob_lidar[idx] / total(obs_prob_lidar)

;-----------------------------------------------------------------------------------------
; smooth
;-----------------------------------------------------------------------------------------
;goto, skip_smooth
smooth_width = 3

;--- new smooth
cond_ratio_clear = smooth_nan(cond_ratio_clear,smooth_width,missing)
cond_ratio_water = smooth_nan(cond_ratio_water,smooth_width,missing)
cond_ratio_ice = smooth_nan(cond_ratio_ice,smooth_width,missing)
post_prob_clear = smooth_nan(post_prob_clear,smooth_width,missing)
post_prob_water = smooth_nan(post_prob_water,smooth_width,missing)
post_prob_ice = smooth_nan(post_prob_ice,smooth_width,missing)

skip_smooth:

;--- normalize obs_prob
idx = where(obs_prob_lidar gt 0,cc)
if (cc gt 0) then obs_prob_lidar[idx] = obs_prob_lidar[idx] / total(obs_prob_lidar)
;print, 'sum of obs_prob_lidar after fill = ', total(obs_prob_lidar)

;--- ensure positive ratios




skip_post_process:

;------------------------------------------------------------------------------------
; apply table
;------------------------------------------------------------------------------------
n = n_elements(phase_true_sfc)
phase_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)])])

;   max_prob = max([post_prob_clear[i,j], post_prob_water[i,j] , post_prob_ice[i,j]])
;   if (post_prob_clear[i,j] eq max_prob) then phase_retrieved[k] = 0
;   if (post_prob_water[i,j] eq max_prob) then phase_retrieved[k] = 1
;   if (post_prob_ice[i,j] eq max_prob) then phase_retrieved[k] = 2

    phase_retrieved[k] = 0
    if (post_prob_clear[i,j] lt 0.50) then begin
         phase_retrieved[k] = 1
         if (post_prob_ice[i,j] gt post_prob_water[i,j]) then begin
            phase_retrieved[k] = 2
         endif
    endif

endfor

idx = where(phase_true_sfc ge 0 and phase_retrieved ge 0, count_all)
idx = where(phase_true_sfc eq 0 and phase_retrieved eq 0, count_agree_clear)
idx = where(phase_true_sfc gt 0 and phase_retrieved gt 0, count_agree_cloud)
idx = where(phase_true_sfc eq 1 and phase_retrieved eq 1, count_agree_water)
idx = where(phase_true_sfc eq 2 and phase_retrieved eq 2, count_agree_ice)

if (count_all gt 0) then begin
  pod_clear = float(count_agree_cloud + count_agree_clear) / count_all
endif else begin
  pod_clear = missing
endelse

if (count_agree_cloud gt 0) then begin
  pod_phase = float(count_agree_water + count_agree_ice) / count_agree_cloud
endif else begin
  pod_phase = missing
endelse

print, 'Results for Surface = ', pod_clear, pod_phase

;------------------------------------------------------------------------------------
; make images
;------------------------------------------------------------------------------------
if (~keyword_set(noplot)) then begin

window, 0, xsize = 600, ysize = 600
load_color_table, 33, 0, 1.0, 249
erase, color = 0

view, post_prob_clear, 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 = 'CLEAR'

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, ['CLEAR'], color = 3, /top, /right, charsize = 2, /clear
legend_local, [string(isfc+1,format="(I02)")], color = 3, /top, /left, charsize = 2, /clear
saveimage, classifier_name+'_post_prob_clear_sfc'+string(isfc,format="(I01)")+'.png',/png

window, 1, xsize = 600, ysize = 600
load_color_table, 33, 0, 1.0, 249
erase, color = 0

view, post_prob_water, 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 = 'WATER'

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, ['WATER'], color = 3, /top, /right, charsize = 2, /clear
legend_local, [string(isfc+1,format="(I02)")], color = 3, /top, /left, charsize = 2, /clear
saveimage, classifier_name+'_post_prob_water_sfc'+string(isfc,format="(I01)")+'.png',/png

window, 2, xsize = 600, ysize = 600
load_color_table, 33, 0, 1.0, 249
erase, color = 0

view, post_prob_ice, 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 = 'ICE'

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, ['ICE'], color = 3, /top, /right, charsize = 2, /clear
legend_local, [string(isfc+1,format="(I02)")], color = 3, /top, /left, charsize = 2, /clear
saveimage, classifier_name+'_post_prob_ice_sfc'+string(isfc,format="(I01)")+'.png',/png

window, 3, xsize = 600, ysize = 600
load_color_table, 33, 0, 1.0, 249
erase, color = 0

view, obs_prob_lidar, 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=max(obs_prob_lidar), /mask_missing, gam_fac=0.1

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, ['OBS_PROB'], color = 3, /top, /right, charsize = 2, /clear
legend_local, [string(isfc+1,format="(I02)")], color = 3, /top, /left, charsize = 2, /clear
saveimage, classifier_name+'_obs_prob_lidar_sfc'+string(isfc,format="(I01)")+'.png',/png

endif

;--------------------------------------------
;
;--------------------------------------------
obs_table[*,*,isfc] = obs_prob_lidar
ice_table[*,*,isfc] = cond_ratio_ice
water_table[*,*,isfc] = cond_ratio_water
clear_table[*,*,isfc] = cond_ratio_clear

cycle_sfc:

endfor  ;sfc type loop

;-------------------------------------------------------
; constrain
;-------------------------------------------------------
idx = where(obs_table lt 0.0 and obs_table ne missing,cc) & if (cc gt 0) then obs_table[idx] = 0.00
idx = where(ice_table lt 0.0 and ice_table ne missing,cc) & if (cc gt 0) then ice_table[idx] = 0.00
idx = where(water_table lt 0.0 and water_table ne missing,cc) & if (cc gt 0) then water_table[idx] = 0.00
idx = where(clear_table lt 0.0 and clear_table ne missing,cc) & if (cc gt 0) then clear_table[idx] = 0.00


;save_filename = sensor_name + "_" + day_string+"_"+classifier_name+".sav"

;save, file = save_filename, nc_header, nbins_x, nsfc, x_min, x_bin, $
;               lat_min, lat_max, $
;               solzen_min, solzen_max, $
;               solglintzen_min, solglintzen_max, $
;               zsfc_min, zsfc_max, $
;               snow_class_min, snow_class_max

  ;-------------------------------------------------
  ; make a structure to hold data and attributes
  ;-------------------------------------------------

  ;-- assign attributes
  @attributes_define.pro
  attributes.classifier_name = classifier_name
  attributes.length_classifier_name = strlen(classifier_name)
  attributes.description = nc_header
  attributes.sensor = sensor_name
  attributes.training_source = training_source
  attributes.timestamp = timestamp()
  attributes.cod_clear_ice_thresh = cod_clear_ice_thresh
  attributes.cod_clear_water_thresh = cod_clear_water_thresh
  attributes.rank = nb_rank
  attributes.nsfc = nsfc
  attributes.nchan_used = nchan_used
  attributes.x_name = xtitle
  attributes.nbins_x = nbins_x
  attributes.x_min = x_min
  attributes.x_bin = x_bin
  attributes.y_name = ytitle
  attributes.nbins_y = nbins_y
  attributes.y_min = y_min
  attributes.y_bin = y_bin
  attributes.z_name = ztitle
  attributes.nbins_z = nbins_z
  attributes.z_min = z_min
  attributes.z_bin = z_bin
  attributes.lat_min = lat_min
  attributes.lat_max = lat_max
  attributes.zen_min = zen_min
  attributes.zen_max = zen_max
  attributes.solzen_min = solzen_min
  attributes.solzen_max = solzen_max
  attributes.solglintzen_min = solglintzen_min
  attributes.solglintzen_max = solglintzen_max
  attributes.solscatang_min = solscatang_min
  attributes.solscatang_max = solscatang_max
  attributes.lunzen_min = lunzen_min
  attributes.lunzen_max = lunzen_max
  attributes.lunglintzen_min = lunglintzen_min
  attributes.lunglintzen_max = lunglintzen_max
  attributes.lunscatang_min = lunscatang_min
  attributes.lunscatang_max = lunscatang_max
  attributes.tpw_min = slant_tpw_min
  attributes.tpw_max = slant_tpw_max
  attributes.tsfc_min = tsfc_min
  attributes.tsfc_max = tsfc_max
  attributes.zsfc_min = zsfc_min
  attributes.zsfc_max = zsfc_max
  attributes.zsfc_std_min = zsfc_std_min
  attributes.zsfc_std_max = zsfc_std_max
  attributes.snow_class_min = snow_class_min
  attributes.snow_class_max = snow_class_max
  attributes.solglint_mask_min = solglint_mask_min
  attributes.solglint_mask_max = solglint_mask_max
  attributes.lunglint_mask_min = lunglint_mask_min
  attributes.lunglint_mask_max = lunglint_mask_max
  attributes.coast_mask_min = coast_mask_min
  attributes.coast_mask_max = coast_mask_max
  attributes.city_mask_min = city_mask_min
  attributes.city_mask_max = city_mask_max
  attributes.moon_illum_frac_min = moon_illum_frac_min
  attributes.moon_illum_frac_max = moon_illum_frac_max
  attributes.rut_solzen_thresh = rut_solzen_thresh

  ;--- assign tables
  @lut_2d_define.pro
  lut.clear = clear_table
  lut.ice = ice_table
  lut.water = water_table
  lut.obs = obs_table
  lut.observation_count = observation_count
  lut.cloud_fraction = cloud_fraction
  lut.ice_fraction  = ice_fraction
  lut.water_fraction  = water_fraction
  lut.on_flag = sfc_type_on_flag
  channels_used.wvl = channel_wvl_used

  ;-----------------------------------------
  ; start test
  ;-----------------------------------------
  idx = where(lut.clear lt 0.0 and lut.clear ne missing,cc)
; print, 'number of non-missing and negative clear values = ', cc

  idx = where(lut.ice lt 0.0 and lut.ice ne missing,cc)
; print, 'number of non-missing and negative ice values = ', cc

  idx = where(lut.water lt 0.0 and lut.water ne missing,cc)
; print, 'number of non-missing and negative water values = ', cc



  ;-----------------------------------------
  ; end test
  ;-----------------------------------------


;----------------------------------------------
; write to netcdf file
;----------------------------------------------
nc_filename = path_dir + sensor_name + "_" + day_string+"_"+classifier_name+".nc"

sd_id = ncdf_create(nc_filename,/CLOBBER,/NETCDF4_FORMAT)

write_2d_table, sd_id, lut, attributes, channels_used, success_flag, /class_att

ncdf_close, sd_id

print, 'write success = ', success_flag

end