Commit 405f6836 authored by Greifeneder Felix's avatar Greifeneder Felix

- added option to specify time-series names (for plot output)

- added error handling in case of an empty time-series
- application of snow mask to the final time-series for better performance
parent ce56fb8c
This diff is collapsed.
......@@ -179,7 +179,7 @@ class GEE_pt(object):
def applysnowmask(image):
tmp = ee.Image(image)
sdiff = tmp.select('VH').subtract(snowref)
wetsnowmap = sdiff.lte(-2.6).focal_mode(100, 'square', 'meters', 3)
wetsnowmap = sdiff.lte(-2.6)#.focal_mode(100, 'square', 'meters', 3)
return (tmp.updateMask(wetsnowmap.eq(0)))
......@@ -249,11 +249,11 @@ class GEE_pt(object):
.combine(gee_s1_track_fltd.select('angle')) \
.combine(gee_s1_track_fltd.select('lia'))
if masksnow == True:
# apply wet snow mask
gee_s1_lin = gee_s1_track_fltd.select('VH').map(tolin)
snowref = ee.Image(10).multiply(gee_s1_lin.reduce(ee.Reducer.intervalMean(5, 100)).log10())
gee_s1_track_fltd = gee_s1_track_fltd.map(applysnowmask)
# if masksnow == True:
# # apply wet snow mask
# gee_s1_lin = gee_s1_track_fltd.select('VH').map(tolin)
# snowref = ee.Image(10).multiply(gee_s1_lin.reduce(ee.Reducer.intervalMean(5, 100)).log10())
# gee_s1_track_fltd = gee_s1_track_fltd.map(applysnowmask)
# create average for buffer area - i.e. compute time-series
gee_s1_mapped = gee_s1_track_fltd.map(createAvg)
......@@ -262,12 +262,18 @@ class GEE_pt(object):
vv_sig0 = 10 * np.log10(
np.array([x['properties']['result']['VV'] for x in tmp['features']], dtype=np.float))
ge_dates = np.array([dt.datetime.strptime(x['id'][17:32], '%Y%m%dT%H%M%S') for x in tmp['features']])
if dual_pol == True:
# get vh
vh_sig0 = 10 * np.log10(
np.array([x['properties']['result']['VH'] for x in tmp['features']], dtype=np.float))
ge_dates = np.array([dt.datetime.strptime(x['id'][17:32], '%Y%m%dT%H%M%S') for x in tmp['features']])
vh_sig0_lin = np.array([x['properties']['result']['VH'] for x in tmp['features']], dtype=np.float)
vh_sig0 = 10 * np.log10(vh_sig0_lin)
if masksnow == True:
snowref = 10 * np.log10(np.mean(vh_sig0_lin[vh_sig0_lin > np.percentile(vh_sig0_lin, 5)]))
snowmask = np.where(vh_sig0 - snowref > -2.6)
vh_sig0 = vh_sig0[snowmask]
vv_sig0 = vv_sig0[snowmask]
ge_dates = ge_dates[snowmask]
if returnLIA == True:
# get lia
......@@ -281,6 +287,9 @@ class GEE_pt(object):
[np.float(x['properties']['count']['VV']) / np.float(x['properties']['totcount']['constant']) for x in
tmp['features']], dtype=np.float)
if masksnow == True:
val_count = val_count[snowmask]
if self.buffer <= 100:
valid = np.where(val_count > 0.01)
else:
......@@ -355,7 +364,7 @@ class GEE_pt(object):
def createAvg(image):
gee_roi = ee.Geometry.Point(self.lon, self.lat).buffer(self.buffer)
reduced_img_data = image.reduceRegion(ee.Reducer.median(), gee_roi, 30)
reduced_img_data = image.reduceRegion(ee.Reducer.median(), gee_roi, 30, tileScale=4)
return ee.Feature(None, {'result': reduced_img_data})
if yearlist == None:
......
......@@ -3,6 +3,7 @@ from GEE_wrappers import GEE_pt
import numpy as np
import os
import matplotlib.pyplot as plt
import pandas as pd
def get_map(minlon, minlat, maxlon, maxlat, outpath,
sampling=100,
......@@ -100,7 +101,8 @@ def get_ts(loc,
footprint=50,
masksnow=True,
calc_anomalies=False,
create_plots=False):
create_plots=False,
names=None):
"""Get S1 soil moisture time-series
Atributes:
......@@ -111,6 +113,7 @@ def get_ts(loc,
masksnow: apply automatic wet snow mask
calc_anomalies: (boolean) calculate anomalies
create_plots: (boolean) generate and save time-series plots to workpath
names: (string or list of strings, optional): list of time-series names
"""
......@@ -119,56 +122,80 @@ def get_ts(loc,
else:
loc = [loc]
if names is not None:
if isinstance(names, list):
print('Name list specified')
else:
names = [names]
sm_ts_out = list()
names_out = list()
cntr = 0
for iloc in loc:
# iterate through the list of points to extract
lon = iloc[0]
lat = iloc[1]
if names is not None:
iname = names[cntr]
else:
iname = None
print('Estimating surface soil moisture for lon: ' + str(lon) + ' lat: ' + str(lat))
# initialize GEE point object
gee_pt_obj = GEE_pt(lon, lat, workpath, buffer=footprint)
sm_ts = gee_pt_obj.extr_SM(tracknr=tracknr, masksnow=masksnow, calc_anomalies=calc_anomalies)
# create plots
if create_plots == True:
if calc_anomalies == False:
# plot s1 soil moisture vs gldas_downscaled
fig, ax = plt.subplots(figsize=(6.5,2.7))
line1, = ax.plot(sm_ts.index, sm_ts,
color='b',
linestyle='-',
marker='+',
label='S1 Soil Moisture',
linewidth=0.2)
ax.set_ylabel('Soil Moisture [%-Vol.]')
plotname = 's1_sm_' + str(lon) + '_' + str(lat) + '.png'
else:
fig, ax = plt.subplots(figsize=(6.5,2.7))
line1, = ax.plot(sm_ts.index, sm_ts['ANOM'].values,
color='r',
linestyle='-',
marker='+',
label='S1 Soil Moisture Anomaly',
linewidth=0.2)
x0 = [sm_ts.index[0], sm_ts.index[-1]]
y0 = [0, 0]
line2, = ax.plot(x0, y0,
color='k',
linestyle='--',
linewidth=0.2)
ax.set_ylabel('Soil Moisture Anomaly [%-Vol.]')
#plt.legend(handles=[line1, line2])
plotname = 's1_sm_anom' + str(lon) + '_' + str(lat) + '.png'
plt.setp(ax.get_xticklabels(), fontsize=6)
plt.savefig(workpath + plotname, dpi=300)
sm_ts_out.append(sm_ts)
gee_pt_obj = None
if isinstance(sm_ts, pd.Series) or isinstance(sm_ts, pd.DataFrame):
# create plots
if create_plots == True:
if calc_anomalies == False:
# plot s1 soil moisture vs gldas_downscaled
fig, ax = plt.subplots(figsize=(6.5,2.7))
line1, = ax.plot(sm_ts.index, sm_ts,
color='b',
linestyle='-',
marker='+',
label='S1 Soil Moisture',
linewidth=0.2)
ax.set_ylabel('Soil Moisture [%-Vol.]')
if iname is None:
plotname = 's1_sm_' + str(lon) + '_' + str(lat) + '.png'
else:
plotname = 's1_sm_' + iname + '.png'
else:
fig, ax = plt.subplots(figsize=(6.5,2.7))
line1, = ax.plot(sm_ts.index, sm_ts['ANOM'].values,
color='r',
linestyle='-',
marker='+',
label='S1 Soil Moisture Anomaly',
linewidth=0.2)
x0 = [sm_ts.index[0], sm_ts.index[-1]]
y0 = [0, 0]
line2, = ax.plot(x0, y0,
color='k',
linestyle='--',
linewidth=0.2)
ax.set_ylabel('Soil Moisture Anomaly [%-Vol.]')
#plt.legend(handles=[line1, line2])
if iname is None:
plotname = 's1_sm_anom_' + str(lon) + '_' + str(lat) + '.png'
else:
plotname = 's1_sm_anom_' + iname + '.png'
plt.setp(ax.get_xticklabels(), fontsize=6)
plt.savefig(workpath + plotname, dpi=300)
sm_ts_out.append(sm_ts)
names_out.append(iname)
return sm_ts_out
gee_pt_obj = None
cntr = cntr + 1
if names is not None:
return (sm_ts_out, names_out)
else:
return sm_ts_out
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment