midraed / water Goto Github PK

Does anyone have an example of how they used the 'cfmask' function in water? I'm not sure at which point in running the functions that I should run the 'cfmask' function.

Thanks

Hello,
I get the error below and need help figuring out what may be wrong.

Error in .local(x, ...) :
RasterLayer has no NA cells (for which to compute a distance)

This is the code used
Energy.Balance <- METRIC.EB(image.DN = image, image.SR = image.SR, plain=TRUE, aoi=aoi, WeatherStation = WeatherStation, ETp.coef = 1.2, MTL=MTLfile, n = 5, sat="L8", alb.coeff = "Olmedo", LST.method = "SW", LAI.method = "metric2010", thermalband=image.DN$thermal.low, Z.om.ws = 0.03, anchors.method = 'flexible')

I am using Landsat 8 images downloaded from ESPA, and plotted correctly as far as I know (see below) Also the weather station data is ok.

plot(image)

plot(image.SR)

Thanks a lot in advance for all your help (I am not that experienced with R in general).

Susana.

If delta rah is big and it > 2, start again

This:

New Landsat Product Identifiers

T​o better identify in which collection a scene is contained, a new Landsat Product ID will be created for Landsat data products. ​The current Scene ID will remain in the visible metadata on EarthExplorer (EE), and in the metadata file (MTL.txt) as a cross-reference for easier metadata searching.​

The proposed Product ID:

LXSS_LLL_PPPRRR_YYYYMMDD_yyyymmdd_CC_TX

L = Landsat (constant)
X = Sensor (“C” = OLI/TIRS Combined, “O”= OLI-only, “T” = TIRS-only, “E” = ETM+, “T” = TM, “M”= MSS)
SS = Satellite (“07” = Landsat 7, “08” = Landsat 8)
LLL = Processing Correction level (“L1T” = precision and terrain, “L1G” = systematic terrain,
“L1S” = systematic)
PPP = WRS path
RRR = WRS row
YYYYMMDD = Acquisition Year (YYYY) / Month (MM) / Day (DD)
yyyymmdd = Processing Year (yyyy) / Month (mm) / Day (dd)
CC = Collection number (“01”,”02”)
TX = Collection Category (”RT” = real-time, “T1” = Tier 1 , or “T2” = Tier 2)

incorporates:

• gap filling by gap mass.
• gap filling by near Landsat scene (date).

I think there is something wrong in the hourlyET function. I've just compared the results with the FAO instruction for hourly data: (i.e. EXAMPLE 19. Determination of ETo with hourly data, http://www.fao.org/docrep/X0490E/x0490e08.htm), and for afternoon example the result is ~0.82 vs. 0.63.

Moreover, I think the documentation needs to be improved to clearly indicate which units are required for radiation (kWh/m2 or MJ/m2), and e.g. that argument height is only used for wind speed corrections (not for all meteorological devices as far as I know).

Anyway, good job and keep on developing!

ALOS DEM is a lot better than SRTM, more updated and with less strange values. Is not available everywhere but we should provide a tool to download and assemble DEM from ALOS data.

Add measurement of uncertainty propagation

I have challenges with the Energy Balance calculation of the water package, following the example code found on: https://cran.rproject.org/web/packages/water/vignettes/Landsat8.html. Using my own data, I read my original landsat image and the surface reflectance image with my area of interests (aoi). But, I get the following error when I use this code.

Energy.Balance <- METRIC.EB(image.DN = image_band, image.SR = image_SR_band,
plain = TRUE, aoi = aoi, n = 5,
WeatherStation = WeatherStation, ETp.coef = 1.2, sat = "L8",
alb.coeff = "Olmedo", LST.method = "SW",
LAI.method = "metric2010", thermalband=image.DN$thermal.low,
Z.om.ws = 0.03, anchors.method = "flexible", MTL = MTLfile)

Error in calcAnchors(image = image.TOAr, Ts = Ts, LAI = LAI, plots = T, :
Automatic selection of hot anchors FAILED

While using raw Landsat 8 data I am getting following error.
Error in names(imagestack) <- stack.names : object 'image_TOA' not found.
I am using water version 0.5.

Hello,
I am trying to analyze the images collected using drone. However, I am not getting any idea on adjusting for the MTL file (used for satellite images), can anyone suggest me something on this.

thanks

La función "read.WSdata" debería leer el número que identifique la posición (número de columna) en la tabla de cada serie de datos, en lugar de exigir definir todas las columnas de forma explícita.

Ejemplo:

read.WSdata(WSdata =("/format2).csv"),
datetime.format = "%Y/%m/%d %H:%M",
columns = c(NA,NA,"datetime", NA,"pp","radiation","RH", NA, NA, NA, NA, "temp",wind)
lat=-33.00513, long= -68.86469, elev=927, height= 2,
MTL=MTLfile)

read.WSdata(WSdata =("/format2).csv"),
datetime.format = "%Y/%m/%d %H:%M",
columns = c(datetime=3,"temp"=12, "RH"=7, "pp"=5,
, "radiation"=6, "wind"=4),
lat=-33.00513, long= -68.86469, elev=927, height= 2,
MTL=MTLfile)

quedando:
columns=c("date", "time", "datetime", "temp", "RH", "pp", "radiation", "wind")
#orden de los argumentos fijo
columns= c(NA,NA, 2,12,7,5,6,4)

Hi,
I'm trying to implement the metric model on another language and I have realized that when estimating the sensible heat the model uses the following algorithm:
Ts.datum <- Ts - (DEM - WeatherStation$elev) * 6.49 / 1000
However, in the literature most authors describe this variable as TsDEM, which has a small difference compared to your implementation:
TsDEM = Ts + 0.0065Δz,
where Δz is the elevation of each pixel minus the elevation of a datum (meters), and corresponds to the term (DEM - WeatherStation$elev) of your equation. However, in the original equation this is added to the surface temperature, but in your algorithm it's subtracted.
Can you please confirm if this was intended to work that way or if it's just a mistake.
Regards

Hello everyone
I have tried to load image.SR for Landsat 8 on data in 2022, this issue is the first time has happened to me. For data in 2021 it did not happen when loading image.SR.

image.SR <- loadImageSR(path=raw_data_folder, aoi=aoi)

Error in .local(.Object, ...) :

Error in .rasterObjectFromFile(x, band = band, objecttype = "RasterLayer", :
Cannot create a RasterLayer object from this file. (file does not exist)

I have tried to solve this issue by using the code below

files <- list.files(path = raw_data_folder

•                 , pattern = "_sr_band+[2-7].tif$", full.names = T)
  

files <- list.files(path = raw_data_folder

• , pattern = "_SR_B+[2-7].tif$", full.names = T)

files <- list.files(path = system.file("March1422", package="water")

•                 , pattern = "_SR_B+[2-7].(TIF|tif)$", full.names = T)
  

stack1 <- list()
for(i in 1:6){

• stack1[i] <- raster(files[i])}
  

Warning messages:
1: In [<-(*tmp*, i, value = raster(files[i])) :
implicit list embedding of S4 objects is deprecated
2: In [<-(*tmp*, i, value = raster(files[i])) :
implicit list embedding of S4 objects is deprecated
3: In [<-(*tmp*, i, value = raster(files[i])) :
implicit list embedding of S4 objects is deprecated
4: In [<-(*tmp*, i, value = raster(files[i])) :
implicit list embedding of S4 objects is deprecated
5: In [<-(*tmp*, i, value = raster(files[i])) :
implicit list embedding of S4 objects is deprecated
6: In [<-(*tmp*, i, value = raster(files[i])) :
implicit list embedding of S4 objects is deprecated

image.SR <- do.call(stack, stack1)
image.SR <- crop(image.SR, aoi)
image.SR <- image.SR / 10000
bandnames <- c("B", "G", "R", "NIR", "SWIR1", "SWIR2")
image.SR <- saveLoadClean(imagestack = image.SR,

•                       stack.names = bandnames, 
  

•                       file = "image.SR", 
  

•                       overwrite=TRUE)
  

Error in saveLoadClean(imagestack = image.SR, stack.names = bandnames, :
could not find function "saveLoadClean"

When i tried to run METRIC.EB, the new issue is happened

Energy.Balance <- METRIC.EB(image.DN = image, image.SR = image.SR,

•                         plain=TRUE, aoi=aoi, n = 5, WeatherStation = WeatherStation, 
  

•                         ETp.coef = 1.2, sat="L8", alb.coeff = "Olmedo", LST.method = "SW", 
  

•                         LAI.method = "metric2010", Z.om.ws = 0.03, MTL = MTLfile)
  

Error in calcAnchors(image = image.TOAr, Ts = Ts, LAI = LAI, plots = T, :
Not enough pixels with the conditions for anchor pixels. I
found 0 cold pixels and 0 hot pixels.
In addition: Warning messages:
1: In calcAnchors(image = image.TOAr, Ts = Ts, LAI = LAI, plots = T, :
anchor method names has changed. Old names (CITRA-MCBx) are
deprecated. Options now include 'best', 'random' and 'flexible'
2: In max(Ts[albedo >= 0.13 & albedo <= 0.15 & NDVI >= 0.1 & NDVI <= :
no non-missing arguments to max; returning -Inf
Can you help me, please?
thank you for your understanding and cooperation
Best regards,

Hi all,
running the function dailyET to estimate 24 h reference evapotranspiration yields NA_real values that gives error in estimating ET.24 as follows.

Error in seq.default(0, ceiling(ETr.daily * 1.5), length.out = 50) :
'to' must be a finite number
In addition: Warning message:
In asMethod(object) :
complete map seems to be NA's -- no selection was made

The weather data file has worked for the entire algorithm and is facing issues in the last command. could anyone please help me urgently on this.

I am using the pacakage of water to emitate ET in a basin. However, I find that the evaporation of water is no data and the other land-cover has the value, I repeat the other images and the result is the same as the former. what should I do if I want get the value about water? Besides, the ESPA products can not be received from USGS. Thanks a lot!

Hey man,
I was wondering if you guys have any plans for adding function for estimation surface reflectance for L8?
I'm using Water V:0.8 and it asks users to download surface reflectance from espa website.

Hello
I am trying to use your water toolbox and metric code for Landsat 8. I am using the following command
Energy.Balance <- METRIC.EB(image.DN = image, image.SR = image.SR,

•                         plain=TRUE, n = 1, WeatherStation = WeatherStation, 
  

•                         ETp.coef = 1.2, sat="L8", alb.coeff = "Olmedo", LST.method = "SW", 
  

•                         LAI.method = "metric2010", Z.om.ws = 0.03,  anchors.method = "flexible", MTL = MTLfile)
  

But it is showing the following error message.

Error in .local(x, ...) :
RasterLayer has no NA cells (for which to compute a distance)

I also tried the advanced option but it also shows the same errors while calculating hot and cold pixels using calcAnchors although I could calculate G and Rn. Can anyone help me??

Hi @midraed, Greetings!
I sometimes face this error while using calcAnchors function:

hot.and.cold <- calcAnchors(image = image.TOAr, Ts = Ts, LAI = LAI, plots = TRUE, albedo = albedo, Z.om = Z.om, n = 1,
anchors.method = "flexible", deltaTemp = 5, verbose = FALSE)
Error in apply(x[i, ], 1, which.min) : dim(X) must have a positive length

I don't exactly know how to solve this issue. Please have a look at this and help me through.
Thank you!

When using the function:

dailyET(Wtst, lat=-18.78062 ,long=-40.98651 , elev=255, ET="ETr")

The function returns NA.

Here is the waterWeatherStation object:https://drive.google.com/drive/folders/0B3f-4ojGAWWXVnZBbDRVMmlNd2s?resourcekey=0-wGwq5U6q_dMojItFbxqQoA&usp=sharing

This may be fixed by changing the last line of the function: return(sum(ET.daily)) by return(sum(ET.daily, na.rm=T))

Thanks

Check water for coherency and grammar

add examples on the most important scripts of water package

Hi there,

I'm trying to follow along with the advanced vignette and am not getting the same result. I've been trying to track down the issue. Not sure I found my main problem but did find a bug.

in calcH the calculation of u200.v is looking for values that do not exist in the provided WeatherStation data. Therefore, it is reverting to the default value of 4 m/s for u200.v. Namely WeatherStation$wind and WeatherStation$height. I've found what I think are the right values WeatherStation$daily$wind_mean and WeatherStation$location$height and can succesfully calculate a u200.v value with them:

Z.om.ws <- 0.03 #from function defaults
u.ws <- WeatherStation$daily$wind_mean * 0.41/log(WeatherStation$location$height/Z.om.ws)
u200.v <- u.ws/0.41 * log(200/Z.om.ws) #returns 6.29

Thanks for all your work on this!

Habria que editar los 3 vignettes (en vignettes/) y editar donde sugieren utilizar ESPA. La idea es sugerir como obtener imagenes de surface reflectance desde earth explorer (lo que antes haciamos con ESPA).

Hi Guilhermo,

I am trying to import my weather station, but it is showing the error:

Error in data.frame(date = unique(WSdata$date), radiation_sum = tapply(WSdata$radiation, :
arguments imply differing number of rows: 1, 0
In addition: Warning message:
In read.WSdata(CSVFILE, date.format = "%d/%m/%Y", lat = -35.42222, :
As tz = "", assuming the weather station time zone is America/Sao_Paulo

Do you know, what is going on? I imported your apple.csv, and works but mine file not.

Thank you,

Error cuando no coinciden temporalmente, debería tirar un warning!

add some optional code to plot on Google Earth based on plotKML

Can anyone please tell when we should use the cloud mask function in the water package? It is not clear at what stage we should use this function.

Test functionality with other dates using L7 and L8

I am using water package to estimate actual ET.However, my study area is larger than your demo data, so I need 2 or 3 landsat images to cover, and more than 6 weather stations located in it. What should I do if I want to get a good result using water package? Thanks a lot!

La función "prepareSRTMdata" debería ser capaz de utilizar la extensión de cualquier archivo espacial. Un ejemplo de esto podría ser el uso del correspondiente al aoi en lugar del de la imagen.

Actual:
prepareSRTMdata(path = getwd(), format = "tif", extent=image)

Ejemplo:

prepareSRTMdata(path = getwd(), format = "tif", extent=aoi)

I'm having problems trying to run the function to read the weather station data. Not even the examples of the package are working in the last actualization v0.6.
When trying to execute these lines of code, which depicts the example of the function in the documentation, there is an error:
csvfile <- system.file("extdata", "apples.csv", package="water")
MTLfile <- system.file("extdata", "L7.MTL.txt", package="water")
WS <- read.WSdata(WSdata = csvfile, date.format = "%d/%m/%Y",lat=-35.42222, long= -71.38639, elev=201, height= 2.2,MTL = MTLfile)

"Error in Ops.data.frame(WS.after, WS.prev) :
‘-’ only defined for equally-sized data frames"

Solve how to use extended characters on names or different encoding.

The calculation of cold anchors in failing showing the following error
Error in calcAnchors(image = image.TOAr, Ts = Ts, LAI = LAI, plots = F, :
Automatic selection of cold anchors FAILED

I have used the flexible method for L8 images. Please suggest on what should be done now.

thanks

I am trying to use the METRIC.EB function (following the Landsat 8 and water package simple procedure vignette), but the MTL file seems to be throwing it off; I get this error: 'Error in readLines(MTL, warn = FALSE) : 'con' is not a connection'.

I read in the MTL file using a 'read.delim' function, and it is properly stored in the R environment, but perhaps this is the wrong way considering it is not a system file within the water package.

Landsat 8 OLI/TIRS and L8 surface reflectance bands and the weather station data were added into the R environment properly ... I can plot all of the bands

Thank you.

Hello

I'm having a issue while reading some weather station data:

csvfile<-"E:/estatistica/projetos/metric/estacoes/ea614.csv"
MTLfile<-"E:/estatistica/projetos/metric/imagens/landsat/LC08_L1TP_216073_20160531_20170324_01_T1/LC08_L1TP_216073_20160531_20170324_01_T1_MTL.txt"

read.WSdata(WSdata =csvfile, height=2, lat=-19.356923 ,long=-40.068680 , elev=38, tz="UTC", date.format = "%Y-%m-%d", time.format = "%h",MTL = MTLfile, columns=c("date","time", "radiation", "wind", NA, "RH", "temp", "pp"))

Error in data.frame(date = unique(WSdata$date), radiation_sum = tapply(WSdata$radiation, :
arguments imply differing number of rows: 1, 0

There is any mistake?

Here is the used data:

https://drive.google.com/open?id=0B3f-4ojGAWWXVnZBbDRVMmlNd2s

Nevermind... The error ocurred because the time format was wrong (%h in place of %H) and the data was with a extra column in front of it from when I used the command write.csv()

It is recommended to create a new error message when the result of this code is all NA:
datetime <- strptime(datetime, format = datetime.format,
tz = tz)

It is recommended to change the way the correct data columns order is informed. It is a little confusing.

Sorry for bothering. Thanks

Add value field with the output, unit and rescaling factor for many functions

I downloaded recent Landsat 8 images that have names of patterns (see below) and the LoadImage() can't load/stack the. Maybe the patterns defined in the function can't detect these names? or I need to rename the files following a certain pattern defined in the f()?

[1] "LC08_L1TP_009052_20200329_20200409_01_T1_2020-03-29_B1.TIF" "LC08_L1TP_009052_20200329_20200409_01_T1_2020-03-29_B10.TIF"
[3] "LC08_L1TP_009052_20200329_20200409_01_T1_2020-03-29_B11.TIF" "LC08_L1TP_009052_20200329_20200409_01_T1_2020-03-29_B2.TIF"
[5] "LC08_L1TP_009052_20200329_20200409_01_T1_2020-03-29_B3.TIF" "LC08_L1TP_009052_20200329_20200409_01_T1_2020-03-29_B4.TIF"
[7] "LC08_L1TP_009052_20200329_20200409_01_T1_2020-03-29_B5.TIF" "LC08_L1TP_009052_20200329_20200409_01_T1_2020-03-29_B6.TIF"
[9] "LC08_L1TP_009052_20200329_20200409_01_T1_2020-03-29_B7.TIF" "LC08_L1TP_009052_20200329_20200409_01_T1_2020-03-29_B8.TIF"
[11] "LC08_L1TP_009052_20200329_20200409_01_T1_2020-03-29_B9.TIF" "LC08_L1TP_009052_20200329_20200409_01_T1_2020-03-29_BQA.TIF"

Just trying to get the "simple" example code to work. The issue is with running:

WeatherStation <- read.WSdata(WSdata = csvfile, date.format = "%d/%m/%Y",
lat=-35.42222, long= -71.38639, elev=201, height= 2.2,
MTL = MTLfile)

I receive the error:

Error in Ops.data.frame(WS.after, WS.prev) :
‘-’ only defined for equally-sized data frames

I am using the latest R version 3.2.4 and RStudio.

Kind regards,

John

I have downloaded all the required grid files in the working directory as suggested by checkSRTMgrids. Now there appears an error when trying to run this function.
prepareSRTMdata()
Error in prepareSRTMdata() :
You need to download SRTM grids and save them to working directory
try: checkSRTMgrids()

I think there is again the issue of data format, as i saw the arguments of the function, it says format ".tif" however, the files are of type "TIF". could you please check and suggest something on this.

thanks

I have challenges with the Energy Balance calculation of the water package, following the example code found on: https://cran.rproject.org/web/packages/water/vignettes/Landsat8.html. Using my own data, I read my original landsat image and the surface reflectance image with my area of interests (aoi). But, I get the following error when I use this code.

Energy.Balance <- METRIC.EB(image.DN = image_band, image.SR = image_SR_band,
plain = TRUE, aoi = aoi, n = 5,
WeatherStation = WeatherStation, ETp.coef = 1.2, sat = "L8",
alb.coeff = "Olmedo", LST.method = "SW",
LAI.method = "metric2010", thermalband=image.DN$thermal.low,
Z.om.ws = 0.03, anchors.method = "flexible", MTL = MTLfile)

Error in calcAnchors(image = image.TOAr, Ts = Ts, LAI = LAI, plots = T, :
Not enough pixels with the conditions for anchor pixels. I
found 0 cold pixels and 0 hot pixels.
In addition: Warning messages:
1: PROJ support is provided by the sf and terra packages among others
2: In calcAnchors(image = image.TOAr, Ts = Ts, LAI = LAI, plots = T, :
anchor method names has changed. Old names (CITRA-MCBx) are
deprecated. Options now include 'best', 'random' and 'flexible'
3: In max(Ts[albedo >= 0.13 & albedo <= 0.15 & NDVI >= 0.1 & NDVI <= :
no non-missing arguments to max; returning -Inf

Change dailyET to use MTL file and estimate dailyET for the same day as the satellite.

The Landsat image in the demo data should be a raw version. In this way, the importing process could be demostrated.

While running the example from the simple procedure vignette, the METRIC.EB function gives incorrect values for the sensible and latent heat, for some pixels.

Here is the call to the function as well as the warning message:

Energy.Balance <- METRIC.EB(image.DN = image.DN, plain=TRUE, 
                              WeatherStation = WeatherStation, 
                              ETp.coef = 1.2, MTL=MTLfile, 
                              sat="L7", thermalband=B6)
# Warning in setValues(r, log(values(x), base = base)) : NaNs produced
#    pixel      X        Y       Ts  LAI type
#1  43077 280560 -3917040 325.4446 0.13  hot
#2 118655 277020 -3921600 311.2627 3.32 cold
# Warning in calcH(anchors = hot.and.cold, Ts = Ts, Z.om = Z.om, WeatherStation = WeatherStation,  :
#   u200 less than threshold value = 0.2942m/s. using u200 = 4m/s

plot(subset(Energy.Balance, c(3,4)))

Here is a zoom on problematic pixels:

out <- crop(Energy.Balance, extent(278000, 280000, -3923000, -3921000))
plot(subset(out, c(3,4)))

Some information:

R version 3.2.3 Patched (2015-12-24 r69814)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 14.04.3 LTS

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
1] water_0.3    rgdal_1.1-3  raster_2.5-2 sp_1.2-1     sos_1.3-8   
[6] brew_1.0-6  

loaded via a namespace (and not attached):
[1] Rcpp_0.12.2     grid_3.2.3      lattice_0.20-33

Hi all, this package is very useful however I am facing some issues during its application. Since I am new with R language, it is a bit difficult for me to fix things.

1. How can I escape the area of interest and load the entire raster image?
2. what are the prerequisites for the digital elevation model to be used during processing the data?
3. I had to move my data to extdata else before that it was showing error "subscript out of bounds"

I would really appreciate if anyone could help me solve this issue.

thank you