mikejohnson51 / climateR

An R 📦 for getting point and gridded climate data by AOI
https://mikejohnson51.github.io/climateR/
MIT License
168 stars 40 forks source link

getTerraClim:"Error in matrix(var, nrow = l, byrow = F) : data is too long #43

Closed DorotheaDeus closed 3 years ago

DorotheaDeus commented 3 years ago

Dear ALL, I am trying to downloading Terraclimate monthly timeseries data using getTerraClim function. This is my code Step TWO

Load packages

library(chirps) library(sf) library(AOI) library(climateR) library(terra) library(raster) library(rasterVis) library(shapefiles) library(sp)

If you get this Error: Requested AOI not in model domain, then deactivate s2 in sf library

library(sf) sf_use_s2(FALSE)

Load babati stn coordinate shapefile as AOI

df<- read.csv('~/Documents/MoW/Year2_2021_2022/Task 2_Capacity_Building/Training_1_RS/Day2/Data/Babati_stn.csv')#read/import station location data show(df)

Set coordinates

coordinates(df)=~ Lon + Lat

Create a shape file and set the coordinate reference system(CRS) as Arc 1960 EPSG 4210

df.shp = df

crs(df.shp)<-'+proj=longlat +datum=WGS84 +no_defs'

crs(df.shp)<-'+proj=longlat +ellps=clrk80 +towgs84=-160,-6,-302,0,0,0,0 +no_defs' #Arc 1960 EPSG 4210

Create as sf object

Create sf extends data.frame-like objects with a simple feature list column

df.shp.sf = as(df.shp, "sf")

plot the point

plot(df.shp.sf)

Downloading Terraclimate monthly timeseries data using getTerraClim function

Babati.ppt=getTerraClim(df.shp.sf, param = 'prcp',startDate = "1918-01-01",endDate = "1918-07-31")

After running I am getting this; "Error in matrix(var, nrow = l, byrow = F) : data is too long In addition: Warning messages: 1: In max(d$year) : no non-missing arguments to max; returning -Inf 2: In min(d$year) : no non-missing arguments to min; returning Inf 3: In min(d$month.index) : no non-missing arguments to min; returning Inf 4: In max(d$month.index) : no non-missing arguments to max; returning -Inf"

Any suggestion will be appreciated

mikejohnson51 commented 3 years ago

@DorotheaDeus, would you be able to share some of the site locations from the CSV? It'll help to know what area of the world this is exploring. Thanks!

DorotheaDeus commented 3 years ago

@DorotheaDeus, would you be able to share some of the site locations from the CSV? It'll help to know what area of the world this is exploring. Thanks!

Thank you for your quick response. I am trying to download just a for single point in a csv file. Sorry I am a newbie I don't know how to attach a file, but the file consist of a single, station located in northern Tanzania, East Africa with the following coordinate;

Babati | Lat | Lon 1 | -4.2078 | 35.7461

I would like to download also for more than one station, unfortunately I don't know how, could you please share with me any a kind of script to download multiple station from a csv or txt file.

I have been also trying to download for the whole country, Tanzania but it is not working, I am getting the same error.

Thank you in advance.

mikejohnson51 commented 3 years ago

Hi, so with the current version of climateR the s2 issues are fixed. Remember that terraclimate only goes back to 1958, so the 1918 date was likely causing your issue. Otherwise, here's a sample to tackle some of your goals:

library(sf)
#> Linking to GEOS 3.9.1, GDAL 3.3.0, PROJ 8.0.1
library(dplyr)
library(raster)
library(climateR)
library(ggplot2)

# Get AOI
AOI = AOI::aoi_get(country = 'Tanzania')

# Get PPT rasters
Babati.ppt = getTerraClim(AOI = AOI::aoi_get(country = 'Tanzania'), 
                        param = 'prcp',
                        startDate = "1958-01-01",
                        endDate = "1958-12-31")
#> Spherical geometry (s2) switched off
#> Spherical geometry (s2) switched on

#Check
plot(Babati.ppt$terraclim_prcp)

# Make some random points (use yours!)
pts = st_sample(AOI, 100) %>% 
  st_as_sf() %>% 
  mutate(ID = 1:n())

# Check
plot(Babati.ppt$terraclim_prcp[[1]])
plot(pts, add = TRUE, pch = 16)

# Extract timeseries
ts = extract_sites(Babati.ppt, pts, "ID")

# Check
ts[[1]] %>% 
  tidyr::pivot_longer(-date) %>% 
  ggplot(aes(x = date, y = value, color = name)) + 
  geom_line() + 
  theme(legend.position = "none")

Created on 2021-08-03 by the reprex package (v2.0.0)

DorotheaDeus commented 3 years ago

Hi, so with the current version of climateR the s2 issues are fixed. Remember that terraclimate only goes back to 1958, so the 1918 date was likely causing your issue. Otherwise, here's a sample to tackle some of your goals:

library(sf)
#> Linking to GEOS 3.9.1, GDAL 3.3.0, PROJ 8.0.1
library(dplyr)
library(raster)
library(climateR)
library(ggplot2)

# Get AOI
AOI = AOI::aoi_get(country = 'Tanzania')

# Get PPT rasters
Babati.ppt = getTerraClim(AOI = AOI::aoi_get(country = 'Tanzania'), 
                        param = 'prcp',
                        startDate = "1958-01-01",
                        endDate = "1958-12-31")
#> Spherical geometry (s2) switched off
#> Spherical geometry (s2) switched on

#Check
plot(Babati.ppt$terraclim_prcp)

# Make some random points (use yours!)
pts = st_sample(AOI, 100) %>% 
  st_as_sf() %>% 
  mutate(ID = 1:n())

# Check
plot(Babati.ppt$terraclim_prcp[[1]])
plot(pts, add = TRUE, pch = 16)

# Extract timeseries
ts = extract_sites(Babati.ppt, pts, "ID")

# Check
ts[[1]] %>% 
  tidyr::pivot_longer(-date) %>% 
  ggplot(aes(x = date, y = value, color = name)) + 
  geom_line() + 
  theme(legend.position = "none")

Created on 2021-08-03 by the reprex package (v2.0.0)

Thank you so much, I appreciate your assistance

DorotheaDeus commented 3 years ago

Hello sorry for bothering, I have this file I was trying to create a simple feature, sf file but I couldn't succeed, STNs.csv then I tried to import it as vector so as to extract the data at multiple sites below Archive.zip but when I try to extract I am getting more than 20 sites, how can I extract correct just the 20 sites I have. How can I associate with the station name instead of station Id. Tried also to plot I could not succeeed. This is how I did,

Screen Shot 2021-08-05 at 15 17 09 Screen Shot 2021-08-05 at 15 17 48

How can I change the color of the raster plotted?

Screen Shot 2021-08-05 at 15 19 06 Screen Shot 2021-08-05 at 15 19 45

I have 20 station points but the ts has more site, I don"t know why? How can I rename them with their station name, these are actual meteorological stations. lastly sorry for bothering I didn't understand the functions used to plot, I tried to search but I couldn't see, please assist

Screen Shot 2021-08-05 at 15 23 29

I would like to put a legend with the station names, how can I change the color of the lines I hope you will find some time to assist. Any idea would be appreciated. Thank you in advance

mikejohnson51 commented 3 years ago

Hi, thanks for documenting your efforts. Please try something like this:

library(ggplot2)
library(dplyr)
library(sf)
#> Linking to GEOS 3.9.1, GDAL 3.3.0, PROJ 8.0.1
library(climateR)

stns = read.csv("/Users/mjohnson/Downloads/STNs.csv") %>% 
  st_as_sf(coords = c('Lon', 'Lat'), crs = 4326)

Babati.ppt = getTerraClim(AOI = stns, 
                          param = 'prcp',
                          startDate = "1958-01-01",
                          endDate = "1958-12-31")
#> Spherical geometry (s2) switched off
#> Spherical geometry (s2) switched on

ext = extract_sites(Babati.ppt, stns, id = 'S_Name')

ext$terraclim_prcp %>% 
  tidyr::pivot_longer(-date) %>% 
  mutate(value = as.numeric(value)) %>% 
  ggplot(aes(x = date, y = value, color = name)) + 
  geom_line() + 
  # Change colors (example)
  scale_color_viridis_d() +
  theme_bw() + 
  theme(legend.position = "bottom")

Created on 2021-08-05 by the reprex package (v2.0.0)

DorotheaDeus commented 3 years ago

Thank you so much