Extract Time Series of a Point ( Lon, Lat) from Netcdf in R

Extract time series of a point ( lon, lat) from netCDF in R

Perhaps use the raster package, this won't work for all NetCDF files but it does for yours:

library(raster)
## brick reads all 22280 layers
r <- brick("obs.nc", varname = "tasmin")
## extract works for all time steps
vals <- extract(r, matrix(c(-120, 52.5), ncol = 2))

dim(vals)
## [1] 1 22280

Note that gives a 1-row, many column matrix because I only gave a single point to extract().

(The extraction is simple with direct copy from the nearest cell, use method = "bilinear" to do interpolation). See ?extract for other options.

Extract time series from netCDF in R

extraction of points data can be easily accomplished by creating a SpatialPoints object containing the point from which you want to extract data, followed by an extract operation.
Concerining the other topics: The "X"s are added because column names can not start with numerals, so a character is added. The horizontal ordering can be easily changed after extraction with some transposing

This, for example, should work (It solves also the "X"s problem and changes the format to "column like"):

library(raster)
library(stringr)
library(lubridate)
library(tidyverse)

b <- brick('/home/lb/Temp/buttami/PRECIPITATION_1998.nc')
lon = c(-91.875,-91.625) # Array of x coordinates
lat <- c(13.875, 14.125) # Array of y coordinates
points <- SpatialPoints(cbind(lat,lon)), # Build a spPoints object

# Etract and tidy
points_data <- b %>%
raster::extract(points, df = T) %>%
gather(date, value, -ID) %>%
spread(ID, value) %>% # Can be skipped if you want a "long" table
mutate(date = ymd(str_sub(names(b),2))) %>%
as_tibble()

points_data

# A tibble: 365 × 3
date `1` `2`
<date> <dbl> <dbl>
1 1998-01-01 0 0
2 1998-01-02 0 0
3 1998-01-03 0 0
4 1998-01-04 0 0
5 1998-01-05 0 0
6 1998-01-06 0 0
7 1998-01-07 0 0
8 1998-01-08 0 0
9 1998-01-09 0 0
10 1998-01-10 0 0
# ... with 355 more rows

plot(points_data$date,points_data$`1`)

Use R to extract time series from netcdf data

This might do it, where "my.variable" is the name of the variable you are interested in:

library(survival)
library(RNetCDF)
library(ncdf)
library(date)

setwd("c:/users/mmiller21/Netcdf")

my.data <- open.nc("my.netCDF.file.nc");

my.time <- var.get.nc(my.data, "time")

n.latitudes <- length(var.get.nc(my.data, "lat"))
n.longitudes <- length(var.get.nc(my.data, "lon"))

n.dates <- trunc(length(my.time))
n.dates

my.object <- var.get.nc(my.data, "my.variable")

my.array <- array(my.object, dim = c(n.latitudes, n.longitudes, n.dates))
my.array[,,1:5]
my.vector <- my.array[1, 2, 1:n.dates] # first latitude, second longitude
head(my.vector)

baseDate <- att.get.nc(my.data, "NC_GLOBAL", "base_date")
bdInt <- as.integer(baseDate[1])

year <- date.mdy(seq(mdy.date(1, 1, bdInt), by = 1,
length.out = length(my.vector)))$year

head(year)

Extracting values for specific lat long from netcdf

You can use the following code for data extraction from .nc files

dat <- structure(list(locatioID = paste0('ID', 1:16), lon = c(73.73, 86, 73.45, 86.41, 85.36, 81.95, 82.57, 75.66, 82.03, 
81.73, 85.66, 85.31, 81.03, 81.70, 87.03, 73.38),
lat = c(24.59, 20.08, 22.61, 23.33, 23.99, 19.09, 18.85, 15.25, 26.78,
16.63, 25.98, 23.28, 24.5, 21.23, 25.08, 21.11)),
row.names = c(1L, 3L, 5L, 8L, 11L, 14L, 17L, 18L, 19L, 21L,
23L, 26L, 29L, 32L, 33L, 35L), class = "data.frame")

temp <- brick("chirps-v2.0.1981.days_p05.nc")

xy <- dat[,2:3] #Column 1 is longitude and column 2 is latitude
xy
spts <- SpatialPoints(xy, proj4string=CRS("+proj=longlat +datum=WGS84"))
#Extract data by spatial point
temp2 <- extract(temp, spts)
temp3 <- t(temp2) #transpose raster object
colnames(temp3) <- dat[,1] #It would be better if you have the location names corresponding to the points
head(temp3)
write.csv(temp3, "Rainfall.csv")

Extracting SST time series at multiple lat, lon locations using CDO

So here is my attempt to answer the question as it was stated in the comments (i.e. you wanted an index which was the midshelf locations averaged and then subtracting the same latitude SST sampled at Longitude=9E). I assume the locations are stored pair-wise in a text file called "locations.txt" as in your question above. The loop part of the answer is from one of this question's solutions.

# first loop over the pairs of indices in the text files.
while read -r -a fields; do
for ((i=0; i < ${#fields[@]}; i += 2)); do
# precise lon/lat for sampled mid-shelf
cdo remapnn,"lon=${fields[i]}/lat=${fields[i+1]}" in.nc pt_${i}.nc
# same lat but lon=9E (change as wanted)
cdo remapnn,"lon=9/lat=${fields[i+1]}" in.nc 9E_${i}.nc
done
done < example.txt

# now take the ensemble average over the points.
cdo ensmean pt_*.nc mid_shelf_sst.nc
cdo ensmean 9E_*.nc mid_shelf_9E.nc

# and calculate the index
cdo sub mid_shelf_sst.nc mid_shelf_9E.nc SST_index.nc


Related Topics



Leave a reply



Submit