-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathreading_extracting_data_nc_files.r
85 lines (66 loc) · 1.56 KB
/
reading_extracting_data_nc_files.r
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
install.packages("ncdf4")
library(ncdf4)
install.packages("plotrix")
library(plotrix)
# read in Chl-a data from copernicus database
fn = "C:/Users/tryst/Desktop/PhD/Chl_a_Baltic/CHLa_SW_BalticSea_2017_07_01-2017_08_31_Copernicusdatabase.nc"
nc = nc_open(fn)
print(nc)
attributes(nc$var)
# format 3d data matrix
# 3D data matrix is longitude:latitude:time
lon = ncvar_get(nc, "lon")
nlon = dim(lon)
head(lon)
lat= ncvar_get(nc, "lat", verbose = FALSE)
nlat = dim(lat)
head(lat)
print(lat)
print(lon)
print(c(nlon, nlat))
t = ncvar_get(nc, "time")
print(t)
tunits = ncatt_get(nc, "time", "units")
nt = dim(t)
nt
tunits
# extract required data
CHL_array = ncvar_get(nc, "CHL")
dlname = ncatt_get(nc, "CHL", "long_name")
dlname <- ncatt_get(nc,"CHL","long_name")
dunits <- ncatt_get(nc,"CHL","units")
fillvalue <- ncatt_get(nc,"CHL","_FillValue")
dim(CHL_array)
print(CHL_array)
dim(CHL_array)
# etract data for Kiel based on coordinates and time/date
kiel = CHL_array[54:58, 127:131, 1:62]
x = sum(!is.na(kiel))
x
se1 = sd(kiel, na.rm = TRUE)/sqrt(x)
mean(kiel, na.rm = TRUE)
sd(kiel, na.rm = TRUE)
se1
# Ahrenshoop data
ahp = CHL_array[173:177, 131:135, 1:62]
y = sum(!is.na(ahp))
y
se2 = sd(ahp, na.rm = TRUE)/sqrt(y)
mean(ahp, na.rm = TRUE)
sd(ahp, na.rm = TRUE)
se2
# Usedom data
use = CHL_array[272:276, 160:164, 1:62]
z = sum(!is.na(ahp))
z
se3 = sd(ahp, na.rm = TRUE)/sqrt(z)
mean(use, na.rm = TRUE)
sd(use, na.rm = TRUE)
se3
print(dat)
dim(dat)
dat[47, 20, 10]
mean(dat, na.rm = TRUE)
dat2 = apply(dat, c(1, 2), mean, na.rm = TRUE)
print(dat2)
mean(dat2, na.rm = TRUE)