[R-sig-Geo] Landsat soil moisture index
Allen McIlwee
amcilwee at westnet.com.au
Sat Nov 9 05:43:55 CET 2013
Hi Sarah
I'm a student at Flinders University in South Australia looking at ways to
support the adaptive land management of grazing lands in the north-west
corner of the State. I am interested in how we can use various
multi-temporal image analysis techniques to enhance our understanding of
ecosystem dynamics and landscape processes.
I have over 100 cloud free Landsat images for one scene I have a field trip
out to next week. There are 5-6 images per year on average dating back to
1986. Each of the images has been processed by the ESPA who use the Landsat
Ecosystem Disturbance Adaptive Processing System (LEDAPS) to calculate
at-surface reflectance. I think the level of processing for these products
is 2A? The thermal band is given as brightness temperature.
I have a subset of the images I am working with at the moment which
represents the single driest image in every year. I am working with a bunch
of indices, all calculated in ERDAS imagine that can act as proxies for the
dynamic processes I am trying to understand. In terms of landscape drying
out and rehydration processes, I want to compare a number of indices
including the Tassel Cap Wetness Index, Normalised Difference Moisture
Index, Normalised Soil moisture Index, and the soil moisture index mentioned
below.
I can import the surface temperature & ndvi images into R by:
rm(list=ls())
library(raster)
library(sp)
library(rgdal)
setwd("D:\\SMI\\")
ndvi <- raster("a_1987_298_l5_ndvi3.img", format="HFA", datatype='FLT4S',
NAflag=0)
projection(ndvi) <- "+proj=utm +zone=53 +south +ellps=WGS84 +datum=WGS84
+units=m +no_defs"
ndvi
class : RasterLayer
dimensions : 6450, 6874, 44337300 (nrow, ncol, ncell)
resolution : 30, 30 (x, y)
extent : 217875, 424095, 6867295, 7060795 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=utm +zone=53 +south +ellps=WGS84 +datum=WGS84 +units=m
+no_defs +towgs84=0,0,0
data source : D:\SMI\a_1987_298_l5_ndvi3.img
names : a_1987_298_l5_ndvi3
values : -45.92834, 49.92658 (min, max)
band6 <- raster("a_1987_298_L5_band6.img", format="HFA", datatype='FLT4S',
NAflag=0)
projection(band6) <- "+proj=utm +zone=53 +south +ellps=WGS84 +datum=WGS84
+units=m +no_defs"
band6
class : RasterLayer
dimensions : 6450, 6874, 44337300 (nrow, ncol, ncell)
resolution : 30, 30 (x, y)
extent : 217875, 424095, 6867295, 7060795 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=utm +zone=53 +south +ellps=WGS84 +datum=WGS84 +units=m
+no_defs +towgs84=0,0,0
data source : D:\SMI\a_1987_298_L5_band6.img
names : a_1987_298_L5_band6
values : 2497, 4259 (min, max)
The frequency distribution of unique NDVI values that I want to calculate a
min and max brightness temperature is
freq(ndvi)
value count
[1,] -46 1
[2,] -33 1
[3,] -31 1
[4,] -28 1
[5,] -25 1
[6,] -23 3
[7,] -21 2
[8,] -19 1
[9,] -18 4
[10,] -15 3
[11,] -13 1
[12,] -12 1
[13,] -11 2
[14,] -10 2
[15,] -9 3
[16,] -8 3
[17,] -7 4
[18,] -6 4
[19,] -5 6
[20,] -3 4
[21,] -2 8
[22,] -1 12
[23,] 0 11386427
[24,] 1 28
[25,] 2 23
[26,] 3 39
[27,] 4 113
[28,] 5 539
[29,] 6 2663
[30,] 7 8683
[31,] 8 25549
[32,] 9 81499
[33,] 10 250318
[34,] 11 1202842
[35,] 12 3572881
[36,] 13 6462049
[37,] 14 7941048
[38,] 15 6749783
[39,] 16 3972814
[40,] 17 1729160
[41,] 18 612375
[42,] 19 204610
[43,] 20 71817
[44,] 21 27397
[45,] 22 14043
[46,] 23 7939
[47,] 24 4418
[48,] 25 2790
[49,] 26 1805
[50,] 27 1158
[51,] 28 773
[52,] 29 487
[53,] 30 339
[54,] 31 206
[55,] 32 163
[56,] 33 109
[57,] 34 74
[58,] 35 74
[59,] 36 38
[60,] 37 42
[61,] 38 30
[62,] 39 21
[63,] 40 19
[64,] 41 19
[65,] 42 3
[66,] 43 6
[67,] 44 4
[68,] 45 4
[69,] 46 3
[70,] 47 3
[71,] 48 2
[72,] 49 2
[73,] 50 1
Not sure where to go from here
Thanks
Allen
-----Original Message-----
From: Sarah Goslee [mailto:sarah.goslee at gmail.com]
Sent: 09 November 2013 08:00
To: Allen McIlwee
Cc: r-sig-geo at r-project.org
Subject: Re: [R-sig-Geo] Landsat soil moisture index
Hi,
Well, you need to start at the beginning.
Where did you get your Landsat images?
What format are they in?
Has someone already calculated NDVI and surface temperature for you, or do
you need to do that? If so, what level of processing have the images already
received?
How have you imported the images into R?
Do you already know how to perform a linear regression in R?
It's not a difficult task programmatically, but there are a number of steps
that you need to consider, and things you must understand before you can get
good results. Remote sensing isn't necessarily an easy field.
The landsat package and associated references may help you get started.
Sarah
On Fri, Nov 8, 2013 at 8:52 AM, Allen McIlwee <amcilwee at westnet.com.au>
wrote:
> Dear list
>
>
>
> WANG et al. 2009 "Satellite remote sensing applications for surface
> soil moisture monitoring" derived a Landsat soil moisture index:
>
>
>
> SMI = (Tsmax - Ts) / (Tsmax - Tsmin)
>
>
>
> where Tsmax, Tsmin are the maximum and minimum surface temperatures
> for a given NDVI value, and Ts is the surface temperature at a given
> pixel for a given NDVI. Values near 1 (low surface temp & low
> greenness) have the highest estimated soil moisture, while values near
> 0 (high surface temp &
> greenness) have the lowest.
>
>
>
> Could anyone please advise me how I can create an image matrix that
> will extract the min and max surface temperatures for each unique NDVI
> value. I then need to apply a linear regression to the max temp vs
> NDVI & min temp vs NDVI scatter plots, and use the corresponding
> slopes and intercepts to calculate Tsmax and Tsmin
>
>
>
> Tsmax = (a1 x NDVI) + b1
>
> TSmin = (a2 x NDVI) + b2
>
>
>
> Hoping this is not such a difficult task for someone mathematically
> inclined
>
>
>
> Thanks
>
> Allen
>
>
--
Sarah Goslee
http://www.functionaldiversity.org
More information about the R-sig-Geo
mailing list