Summary

In this notebook we will setup a TSEB model to run a single image using local meteorological data.

TSEB setup

The main input data consists of any GDAL compatible raster image containing the radiometric temperature(s). Depending on the TSEB model to run the image should contain:

Priestley-Taylor TSEB (TSEB-PT): 1st band with the radiometric surface temperature (Kelvin)
Dual-Time Difference TSEB (DTD): 1st band with the radiometric surface temperature around noon. 2nd band with radiometric surface temperature around sunrise (Kelvin)
Composite temperatures TSEB( TSEB-2T): 1st band with canopy temperature. 2nd band with soil temperature(Kelvin)

You can add additional images for leaf area index, view zenith angle (degrees), fractional cover, canopy height (m), canopy width-to-height ratio, green fraction of vegetation, and a processing mask. For all these inputs, you can also type a constant value to be applied to the whole scene (i.e. use this constant value to run TSEB in all the pixels within the processing mask).

For more information about the inputs and outpus of pyTSEB click here

Select the input and output files and fill the information in all tabs

You can press the Load Configuration File button to load a configuration text file which will upload all its information to the corresponding cells.

You can also press the Save Configuration File button to save the configuration in a text file which can be used in future runs.



In [ ]:

    
from pyTSEB.TSEBIPythonInterface import TSEBIPythonInterface # import the PyTSEB class object in the pyTSEB.py module
setup = TSEBIPythonInterface() # Create the setup instance from the PyTSEB class object
setup.local_image_widget() # Load the TSEB configuration Widget

Run TSEB

Once TSEB is configured we parsed all the information in the widgets and run the selected model.

After the processing is done two GeoTiff files will be saved:

< Main Output File > whose name is specified in the cell Output File will contain the bulk estimated fluxes:
1. Net radiation (W m-2)
2. Sensible heat flux (W m-2)
3. Latent heat flux (W m-2)
4. Soil heat flux (W m-2)

< Ancillary Output File > with the same name as the main input file but with a suffix _ancillary added, will contain ancillary information from TSEB:
1. Net shortwave radiation (W m-2)
2. Net longwave radiation (W m-2)
3. Canopy sensible heat flux (W m-2)
4. Canopy latent heat flux (W m-2)
5. Evapotrasnpiration partitioning (canopy LE/total LE)
6. Canopy temperature (K)
7. Soil temperature (K)
8. Aerodynamic resistance (s m-1)
9. Bulk canopy resistance to heat transport (s m-1)
10. Soil resistance to heat transport (s m-1)
11. Friction velocity (m s-1)
12. Monin-Obukhov lenght (m)
13. Friction velocity (m s-1)
14. Quality Flag (unitless)

Display results

Now we can open the image and display the TSEB outputs:



In [ ]:

    
# Change to have a different colour stretch
high_flux=600 # Maximum flux value in the display
low_flux=0 # Minimum flux value in the display

from bokeh.plotting import *
from bokeh.palettes import RdYlBu11 as colortable
from bokeh.models.mappers import LinearColorMapper
from bokeh.io import output_notebook
import numpy as np
import gdal
from bokeh.resources import INLINE
output_notebook(resources=INLINE)

# Open the file
fid=gdal.Open(setup.outputFile,gdal.GA_ReadOnly)
rows=fid.RasterYSize
cols=fid.RasterXSize
# read each band and store the arrays
H=fid.GetRasterBand(1).ReadAsArray()
LE=fid.GetRasterBand(2).ReadAsArray()
Rn=fid.GetRasterBand(3).ReadAsArray()
G=fid.GetRasterBand(4).ReadAsArray()
del fid
colortable=list(reversed(colortable))
map_LE=LinearColorMapper(palette=colortable,high=high_flux,low=low_flux)

# Setup the figure
s1= figure(title="H",plot_width=cols, plot_height=rows, x_range=[0, cols], y_range=[0, rows])
s1.axis.visible = None
s1.image(image=[np.flipud(H)],x=[0],y=[0],dw=cols,dh=rows,color_mapper=map_LE)
s2= figure(title="LE",plot_width=cols, plot_height=rows, x_range=s1.x_range, y_range=s1.y_range)
s2.axis.visible = None
s2.image(image=[np.flipud(LE)],x=[0],y=[0],dw=[cols],dh=[rows],color_mapper=map_LE)
s3= figure(title="Rn",plot_width=cols, plot_height=rows, x_range=s1.x_range, y_range=s1.y_range)
s3.image(image=[np.flipud(Rn)],x=[0],y=[0],dw=[cols],dh=[rows],color_mapper=map_LE)
s3.axis.visible = None
s4= figure(title="G",plot_width=cols, plot_height=rows, x_range=s1.x_range, y_range=s1.y_range)
s4.image(image=[np.flipud(G)],x=[0],y=[0],dw=[cols],dh=[rows],color_mapper=map_LE)
s4.axis.visible = None
p = gridplot([[s1, s2,s3,s4]], toolbar_location='above')

# Add a colormap legend
y = np.linspace(low_flux,high_flux,len(colortable))
dy = y[1]-y[0]
ramp = figure(tools="", y_range = [0, 1], x_range = [low_flux,high_flux], plot_width = 650, plot_height=100)
ramp.toolbar_location=None
ramp.yaxis.visible = None
ramp.rect(x=y, y=0.5, color=colortable, width=dy, height = 1)

show(p)
show(ramp);