;"$Id:$" ;--------------------------------------------------------------- @read_group_attributes @read_1d_group_lut @read_2d_group_lut @read_3d_group_lut @read_prior @get_prob_1d @get_prob_2d @get_prob_3d pro read_combined_table, nc_filename Missing = -999.0 Prob_Min_Thresh = 0.001 R_Max = 1000.0 ;------------------------------------------------------------------- ; read lookup table ;------------------------------------------------------------------- ;--- open file for output sd_id = ncdf_open(nc_filename) ;--- open file and get the number and names of classifiers ncdf_attget,sd_id, "nclassifiers", nclass, /global ;-- open lut file sds_id = ncdf_varid(sd_id,'classifier_names') ncdf_varget, sd_id, sds_id, classifier_names_tmp ;--- convert byte array to string array classifier_names = string(classifier_names_tmp) ;--- define space lut=ptrarr(nclass, /allocate_heap) attrs=ptrarr(nclass, /allocate_heap) channels=ptrarr(nclass, /allocate_heap) ;--- get group id numbers - each group is a classifier grp_ids = ncdf_groupsinq(sd_id) class_read_flag = make_array(Nclass,/INT,VALUE=0) for iclass = 0, Nclass -1 do begin grp_id = grp_ids[iclass] classifier_name = ncdf_groupname(grp_id) ;--- read attributes read_group_attributes, grp_id, attrs_temp *(attrs[iclass]) = attrs_temp case attrs_temp.rank of 1: read_1d_group_lut, grp_id, lut_temp, attrs_temp, channels_temp, success_flag 2: read_2d_group_lut, grp_id, lut_temp, attrs_temp, channels_temp, success_flag 3: read_3d_group_lut, grp_id, lut_temp, attrs_temp, channels_temp, success_flag else: begin print, 'unknown rank of classifier in read_combined_table' success_flag = 0 end endcase *(lut[iclass]) = lut_temp *(channels[iclass]) = channels_temp class_read_flag[iclass] = success_flag endfor ;--- close file lut file ncdf_close, sd_id stop ;------------------------------------------------------------------- ; read prior table ;------------------------------------------------------------------- Prior_Name = 'prior/nb_cloud_mask_modis_prior.nc' read_prior, Prior_Name, Prior ;--------------------------------------------------------------------------------------------- ; read level2 ;--------------------------------------------------------------------------------------------- L2_Name = 'level2/clavrx_MYD02SSH.A2011078.1310.006.2014079205419.level2.hdf' read_l2, L2_Name, L2 ;--------------------------------------------------------------------------------------------- ; compute prior ;--------------------------------------------------------------------------------------------- compute_prior_from_map, l2.Lon, l2.Lat, l2.Bad_Pixel_Mask, l2.Nx, l2.Ny, Missing, $ prior.lon_min, prior.dLon, prior.Nlon, $ prior.lat_min, prior.dLat, prior.Nlat, $ prior.table, Pixel_Prior ;--------------------------------------------------------------------------------------------- ; run data from ecm ;--------------------------------------------------------------------------------------------- Clear_Cond_Ratio = make_array(Nclass,/FLOAT,Value=Missing) Water_Cond_Ratio = make_array(Nclass,/FLOAT,Value=Missing) Water_Cond_Ratio = make_array(Nclass,/FLOAT,Value=Missing) Obs_Prob = make_array(Nclass,/FLOAT,Value=Missing) Post_Prob_Cloud = make_array(L2.Nx,L2.Ny,/Float,Value=Missing) Post_Prob_Cloud_By_Class = make_array(L2.Nx,L2.Ny,Nclass,/Float,Value=Missing) Post_Prob_Clear = make_array(L2.Nx,L2.Ny,/Float,Value=Missing) Post_Prob_Water = make_array(L2.Nx,L2.Ny,/Float,Value=Missing) Post_Prob_Ice = make_array(L2.Nx,L2.Ny,/Float,Value=Missing) for Elem_Idx = 0, L2.Nx - 1 do begin for Line_Idx = 0, L2.Ny - 1 do begin ;---- reset table probabilities to missing Clear_Cond_Ratio[*] = Missing Water_Cond_Ratio[*] = Missing Ice_Cond_Ratio[*] = Missing Obs_Prob[*] = Missing R_Clear = 1.0 R_Water = 1.0 R_Ice = 1.0 for Class_Idx = 0, Nclasses - 1 do begin ;--- set up input for i = 0, *(attrs[Class_Idx]).rank -1 do begin case i of 0: dim_name = *(attrs[Class_Idx]).x_name 1: dim_name = *(attrs[Class_Idx]).y_name 2: dim_name = *(attrs[Class_Idx]).z_name endcase case dim_name of 'etropo': value = L2.Emiss_Tropo_11_0um_nom[Elem_Idx,Line_Idx] 'topa': value = L2.Topa[Elem_Idx,Line_Idx] 'btd8511': value = L2.btd8511[Elem_Idx,Line_Idx] 'logbt11std': value = L2.logbt11std[Elem_Idx,Line_Idx] endcase case i of 0: x = value 1: y = value 2: z = value endcase endfor ;---- compute for this class if (*(attrs[Class_Idx]).rank eq 1) then begin get_prob_1d, x, l2.solzen[i,j], l2.glintzen[i,j], l2.zsfc[i,j], l2.snow_class[i,j], isfc,$ *(lut[Class_Idx]), *(attrs[Class_Idx]), Prob_Min_Thresh, R_Max, Missing, $ z1, z2, z3, z4 endif if (*(attrs[Class_Idx]).rank eq 2) then begin get_prob_2d, x, y, solzen[i,j], glintzen[i,j], zsfc[i,j], snow_class[i,j], isfc,$ *(lut[itable]), *(attrs[itable]), prob_min_thresh, r_max, missing, $ z1, z2, z3, z4 endif if (*(attrs[Class_Idx]).rank eq 3) then begin get_prob_3d, x, y, z, solzen[i,j], glintzen[i,j], zsfc[i,j], snow_class[i,j], isfc,$ *(lut[itable]), *(attrs[itable]), prob_min_thresh, r_max, missing, $ z1, z2, z3, z4 endif Clear_Cond_Ratio[Class_Idx] = z1 Water_Cond_Ratio[Class_Idx] = z2 Ice_Cond_Ratio[Class_Idx] = z3 Obs_Prob[Class_Idx] = z4 ;--- combine this class into the running product of all classes if (clear_cond_ratio[Class_Idx] ne Missing) then R_Clear = R_Clear * Clear_Cond_Ratio[Class_Idx] if (water_cond_ratio[Class_Idx] ne Missing) then R_Water = R_Water * Water_Cond_Ratio[Class_Idx] if (ice_cond_ratio[Class_Idx] ne Missing) then R_Ice = R_Ice * Ice_Cond_Ratio[Class_Idx] endfor ;-------------------------------------------------------------------------------------- ; compute posterior probs for this pixel ;-------------------------------------------------------------------------------------- ;--- turn ratios into probabilities via post_prob = 1.0 / ( 1.0 + r/prior_yes - r) Post_Prob_Clear[i,j] = 1.0 / (1.0 + R_Clear/Prior_Prob_Clear - R_Clear) Post_Prob_Water[i,j] = 1.0 / (1.0 + R_Water/Prior_Prob_Water - R_Water) Post_Prob_Ice[i,j] = 1.0 / (1.0 + R_Ice/Prior_Prob_Ice - R_Ice) ;--- normalize these probabilties so that they sum to 1 Post_Sum = Post_Prob_clear[i,j] + Post_Prob_water[i,j] + Post_Prob_ice[i,j] Post_Prob_Clear[i,j] = Post_Prob_Clear[i,j] / Post_Sum Post_Prob_Water[i,j] = Post_Prob_Water[i,j] / Post_Sum Post_Prob_Ice[i,j] = Post_Prob_Ice[i,j] / Post_Sum ;--- make a cloud probability from the complement of the clear probability Post_Prob_Cloud[i,j] = 1.0 - Post_Prob_clear[i,j] endfor endfor stop end