Example 1: Use of surface model Oh92 and canopy model SSRT

1. Requirements

• Installation of SenSE

2. Oh92+SSRT retrieval of soil moisture

import numpy as np
#from sense.surface import Dubois95, Oh92
from sense.util import f2lam
from sense.model import RTModel
from sense.soil import Soil
from sense.canopy import OneLayer
import matplotlib.pyplot as plt
import random
from sense.surface import Oh92, Oh04
from scipy.optimize import minimize
import pdb

#### Choose models
#-----------------
canopy = 'turbid_isotropic'
surface = 'Oh92'
models = {'surface' : surface, 'canopy' : canopy}
pol='vv'

# model parameter Oh92
#-----------------
freq = 5.405
clay = 0.0738
sand = 0.2408
bulk = 1.45
theta = np.deg2rad(35)
s = 0.013
sm = np.random.uniform(low=0.05, high=0.35, size=(50,))

# model parameter SSRT
#-----------------
d = 0.55
tau = 0.45
ke = tau/d
omega = 0.175
ks=omega*ke

# run model to produce backscatter
#-----------------
S = Soil(f=freq, s=s, mv=sm, sand=sand, clay=clay, bulk=bulk)
C = OneLayer(ke_h=ke, ke_v=ke, d=d, ks_v=ks, ks_h=ks, canopy=models['canopy'])
RT = RTModel(theta=theta, models=models, surface=S, canopy=C, freq=freq)
RT.sigma0()
back_vv = RT.stot['vv']
back_hv = RT.stot['hv']

# helper function retrieval
#-----------------
def run_model(dic, models):
    # surface
    soil = Soil(mv=dic['mv'], s=dic['s'], clay=dic['clay'], sand=dic['sand'], f=dic['f'], bulk=dic['bulk'])

    # canopy
    can = OneLayer(canopy=dic['canopy'], ke_h=dic['ke'], ke_v=dic['ke'], d=dic['d'], ks_h = dic['omega']*dic['ke'],
                   ks_v = dic['omega']*dic['ke'])

    S = RTModel(surface=soil, canopy=can, models=models, theta=dic['theta'], freq=dic['f'])
    S.sigma0()
    return S.__dict__['stot']['vv'[::-1]], S.__dict__['stot']['vh'[::-1]]

def solve_fun(VALS, var_opt, dic, models):

    for i in range(len(var_opt)):
        dic[var_opt[i]] = VALS[i]

    vv, vh = run_model(dic, models)

    return vv, vh

def fun_opt(VALS, var_opt, dic, models, pol):
    if pol == 'vv':
        return(np.nansum(np.square(solve_fun(VALS, var_opt, dic, models)[0]-dic['vv'])))
    elif pol == 'vh':
        return(np.nansum(np.square(solve_fun(VALS, var_opt, dic, models)[1]-dic['vh'])))
    elif pol == 'vv_vh':
        return(np.nansum(np.square((solve_fun(VALS, var_opt, dic, models)[0]-dic['vv'])/2+(solve_fun(VALS, var_opt, dic, models)[1]-dic['vh'])/2)))

# run soil moisture retrieval
#-----------------

dic = {"mv":0.2, "s":s, "clay":clay, "sand":sand, "f":freq, "bulk":bulk, "canopy":canopy, "d":d,
       "ke":ke, "vv":back_vv, "vh":back_hv, "theta":theta, "omega": omega}

var_opt = ['mv']
guess = [0.2]
bounds = [(0.05,0.35)]

method = 'L-BFGS-B'

sm_retrieved = []

for i,ii in enumerate(back_vv):

    dic = {"mv":0.2, "s":s, "clay":clay, "sand":sand, "f":freq, "bulk":bulk, "canopy":canopy, "d":d,
       "ke":ke, "vv":back_vv[i], "vh":back_hv[i], "theta":theta, "omega": omega}

    res = minimize(fun_opt,guess,args=(var_opt, dic, models, pol),bounds=bounds, method=method)

    fun_opt(res.x, var_opt, dic, models, 'vv')
    sm_retrieved.append(res.x[0])

diff = sm - sm_retrieved
diff_average = np.sum(abs(diff))/len(diff)
print(diff)
print(diff_average)

[ 6.59330364e-06  3.07387396e-05 -1.64056693e-05  1.26983242e-06
  1.59499940e-05  6.26880566e-06  1.51476691e-04  1.82387528e-05
  1.19301459e-06  3.88467200e-05  4.58888086e-06  8.67545295e-06
  2.53841100e-05  1.62021203e-05  4.35594912e-06  1.95046970e-05
  1.85616690e-05  1.74444353e-04  4.85911293e-05 -5.68510333e-06
 -9.97495284e-06  4.40446923e-05  2.19569464e-05 -6.49517331e-06
  3.60644372e-05  2.20953377e-06 -8.19768433e-06 -1.48786469e-06
  1.16335317e-04  1.52661750e-05 -6.06365025e-06 -6.37658096e-06
 -9.41156403e-06  1.37811644e-05  1.08180119e-04  1.25106158e-06
 -1.26556797e-05 -1.28062762e-05  2.77727528e-05 -8.58538885e-06
 -8.22013037e-06  1.86317181e-05 -3.80880104e-06  1.86729943e-06
 -1.27028686e-05  5.38805515e-05  6.44892535e-06  2.94535935e-05
  4.78866799e-06 -1.66380938e-06]
2.446716736412141e-05