Commit 1ee55da6 authored by Greifeneder Felix's avatar Greifeneder Felix

Merge branch 'develop' into 'master'

Develop

See merge request !3
parents b5a48a52 5a3ea15c
This diff is collapsed.
......@@ -14,8 +14,16 @@ from pytesmo.time_series.plotting import plot_clim_anom
from sklearn.linear_model import LinearRegression
from pytesmo.temporal_matching import df_match
# Class to access information from GEE, related to point coordinates TODO include sopport for multiple points
class GEE_pt(object):
"""Class to create an interface with GEE for the extraction of parameter time series
Attributes:
lon: longitude in decimal degrees
lat: latitude in decimal degrees
workdir: path to a directory to save output e.g. time-series plots
buffer: radius of the time-series footprint
"""
def __init__(self, lon, lat, workdir, buffer=20):
ee.Reset()
......@@ -171,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)))
......@@ -241,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)
......@@ -254,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
......@@ -273,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:
......@@ -347,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,8 +3,10 @@ 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,
def get_map(minlon, minlat, maxlon, maxlat,
outpath,
sampling=100,
year=None, month=None, day=None,
tracknr=None,
......@@ -12,6 +14,22 @@ def get_map(minlon, minlat, maxlon, maxlat, outpath,
tempfilter=True,
mask='Globcover',
masksnow=True):
"""Get S1 soil moisture map
Atributes:
minlon, minlat, maxlon, maxlat: (float) extent of soil moisture map
outpath: (string) destination for soil moisture map export
sampling: (integer) map sampling
year, month, day (optional): specify date for extracted map - if not specified, the entire
time-series of maps will be extracted; if specified, the S1 acquisition
closest (in time) will be selected
tracknr (optional): (integer) Use data from a specific Sentinel-1 track only
tempfilter: (boolean) apply multi-temporal speckle filter
mask: (string) select 'Globcover' or 'Corine' (europe only) land-cover classifications for the
masking of the map
masksnow: (boolean) apply snow mask
"""
maskcorine=False
maskglobcover=False
......@@ -94,54 +112,107 @@ def get_map(minlon, minlat, maxlon, maxlat, outpath,
GEE_interface = None
def get_ts(lon, lat,
def get_ts(loc,
workpath,
tracknr=None,
footprint=50,
masksnow=True,
calc_anomalies=False,
create_plots=False):
"""Get S1 soil moisture time-series"""
# 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'
create_plots=False,
names=None):
"""Get S1 soil moisture time-series
Atributes:
loc: (tuple or list of tuples) longitude and latitude in decimal degrees
workpath: destination for output files
tracknr (optional): Use data from a specific Sentinel-1 track only
footprint: time-series footprint
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
"""
if isinstance(loc, list):
print('list')
else:
loc = [loc]
if names is not None:
if isinstance(names, list):
print('Name list specified')
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)
return sm_ts
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)
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)
gee_pt_obj = None
cntr = cntr + 1
if names is not None:
return (sm_ts_out, names_out)
else:
return sm_ts_out
from setuptools import setup
setup(name='pysmm',
version='0.5',
version='0.5.1',
description='Python Sentinel Soil-Moisture Mapping',
long_description='For detailed documentation go-to http://pysmm.readthedocs.io/en/latest/',
classifiers=[
......
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