import numpy as np from numpy.lib.stride_tricks import sliding_window_view _bad_data = -999.0 def spatial_var(rad, threshold): test = sliding_window_view(np.pad(rad, [1, 1], mode='constant'), (3, 3)) - np.expand_dims(rad, (2, 3)) test = np.abs(test.reshape(rad.shape[0], rad.shape[1], 9)) var = np.ones(test.shape) var[test < threshold] = 0 return var.sum(axis=2) def sunglint(viirs_data, threshold, bit, conf): m09 = viirs_data.M02.values m20 = viirs_data.M12.values m31 = viirs_data.M15.values m02 = viirs_data.M07.values m01 = viirs_data.M05.values irclr = np.ones(viirs_data.M02.shape) idx = np.nonzero((m02 > 1.3) & (m01 <= 2.0)) m02[idx] = 1.3 idx = np.nonzero((m02 < -99) | (m02 > 1.3)) m02[idx] = _bad_data irclr[bit < 3] = 0 var = spatial_var(m31, 0.4) reg_var_mean = sliding_window_view(np.pad(m02, [1, 1], mode='constant'), (3, 3)).reshape(m02.shape[0], m02.shape[1], 9).mean(axis=2) reg_std = sliding_window_view(np.pad(m02, [1, 1], mode='constant'), (3, 3)).reshape(m02.shape[0], m02.shape[1], 9).std(axis=2) d37_11 = m20 - m31 idx = np.nonzero((var == 0) & (m09 != _bad_data) & (m20 != _bad_data) & (m31 != _bad_data) & (irclr == 1) & (m02 != _bad_data) & (m09 < threshold['sngm09']) & (d37_11 >= threshold['sg_tbdfl'])) conf[idx] = 0.67 idx = np.nonzero((var == 0) & (m09 != _bad_data) & (m20 != _bad_data) & (m31 != _bad_data) & (irclr == 1) & (m02 != _bad_data) & (m09 < threshold['sngm09']) & (d37_11 >= threshold['sg_tbdfl']) & (reg_var_mean != _bad_data) & (reg_var_mean*reg_std < threshold['sngm02vm'])) conf[idx] = 0.96 return conf def spatial(viirs_data, threshold, scene, conf): m02 = viirs_data.M07.values m31 = viirs_data.M15.values var_m31 = spatial_var(m31, 0.40) var_m02 = spatial_var(m02, 0.0020) idx = np.nonzero((conf > 0.95) & (scene['day'] == 1) & (var_m31 == 0) & (var_m02 == 0)) conf[idx] = 1 idx = np.nonzero((conf > 0.66) & (conf <= 0.95) & (scene['day'] == 1) & (var_m31 == 0) & (var_m02 == 0)) conf[idx] = 0.96 idx = np.nonzero((conf <= 0.66) & (scene['day'] == 1) & (var_m31 == 0)) conf[idx] = 0.67 return conf def land(viirs_data, threshold, scene, conf): m04 = viirs_data.M04.values.ravel() m05 = viirs_data.M08.values.ravel() m20 = viirs_data.M12.values.ravel() m22 = viirs_data.M13.values.ravel() m31 = viirs_data.M15.values.ravel() eco = viirs_data.eco.values.ravel() desert = scene['desert'].ravel() conf = conf.ravel() tbadj = 0 ldsbt11bd = np.array(threshold['Land_Restoral']['ldsbt11bd']) ldsbt11 = np.array(threshold['Land_Restoral']['ldsbt11bd']) irclr = 1 hds11 = np.ones((eco.shape[0], 3)) * (ldsbt11 - tbadj) hds11[eco == 8, :] = ldsbt11bd - tbadj if irclr == 1: conf[m31 > hds11[:, 2]] = 1 conf[(m31 > hds11[:, 1]) & (m31 <= hds11[:, 2])] = 0.96 conf[m31 <= hds11[:, 1]] = 0.5 m5_4_thr = np.full(eco.shape, threshold['Land_Restoral']['ldr5_4_thr']) m5_4_thr[desert == 1] = threshold['Land_Restoral']['ldsr5_4_thr'] m5_4 = m05/m04 md1 = m20 - m22 md2 = m22 - m31 idx = np.nonzero((md1 < threshold['Land_Restoral']['ld20m22']) & (md2 < threshold['Land_Restoral']['ld22m31']) & (m5_4 > m5_4_thr) & (conf <= 0.95)) conf[idx] = 0.96 conf = conf.reshape(viirs_data.M01.shape) return conf def coast(viirs_data, threshold, scene, conf): m01 = viirs_data.M05.values m02 = viirs_data.M07.values coast_ndvi = threshold['Coastal_NDVI_Thresholds']['coast_ndvi'] irclr = 1 if irclr == 1: ndvi = (m02 - m01)/(m01 + m02) idx = np.nonzero((ndvi <= coast_ndvi[0]) | (ndvi >= coast_ndvi[1])) conf[idx] = 1 return conf