Introduction

This paper is meant to demonstrate how to use a simple Kriging interpolation from the RGeostats package applied to the Annual CTD datasets from Norwegian Institute of Marine Research (IMR).

The data set are NetCDF files (*.nc) with following variables: * LATITUDE: Latitude (degrees) of the vessel position. 1D size = {Nb Positions} * LONGITUDE: Longitude (degrees) of the vessel position. 1D size = {Nb Positions} * TIME: Time of the position (number of days since January 1st, 1950). 1D size = {Nb Positions} * DEPH: Depth of each measure (m). 2D size = {Nb Positions, Nb Measures} * TEMP: Temperature (°C). 2D size = {Nb Positions, Nb Measures} * CNDC: Electrical conductivity (S m-1). 2D size = {Nb Positions, Nb Measures} * PSAL: Sea water practical salinity (0.001). 2D size = {Nb Positions, Nb Measures}

Definition of the environment

The next cells have specific contents that the user must choose to run or to skip. Their order is important.

  • Defining the environment for knitr

  • Loading the library Intaros

library(RIntaros)
  • Cleaning the workspace: this paragraph is not systematically performed.
rm(list=ls())
  • Defining if the data set must be read or not from the CSV (flag.read)
flag.read = FALSE

Loading Data

First of all, we setup some environment variables (data file name and bounding box). The flag_file allows the user to store each generated graphic file as a PNG file in the image_name directory, instead of plotting them.

# Setup environment
dir.name   = getwd()
file.name  = "imr_data_0_to_100m.csv"
data.name  = "data"
image.name = "images"
long_lim   = c(-20,60)
lat_lim    = c(52,83)
intaros.save.environ(long_lim = long_lim, lat_lim = lat_lim,
                     flag_file = FALSE,image_name = file.path(dir.name,image.name))

Then we read the CSV file (taking the header line into account) and create the RGeostats Db. Finally we show the contents of the newly created Db (named db0).

if (flag.read || ! exists("db0")) db0 = imr_read_csv(file.path(dir.name,data.name,file.name))
db0
## 
## Data Base Characteristics
## =========================
## 
## Data Base Summary
## -----------------
## File is organized as a set of isolated points
## Space dimension              = 2
## Number of fields             = 10
## Maximum Number of attributes = 10
## Total number of samples      = 5554585
## 
## Variables
## ---------
## Field =   1 - Name    =  rank - Locator =  NA
## Field =   2 - Name    =  Longitude - Locator =  x1
## Field =   3 - Name    =  Latitude - Locator =  x2
## Field =   4 - Name    =  Depth - Locator =  NA
## Field =   5 - Name    =  Temperature - Locator =  NA
## Field =   6 - Name    =  Conductivity - Locator =  NA
## Field =   7 - Name    =  Salinity - Locator =  NA
## Field =   8 - Name    =  Time - Locator =  NA
## Field =   9 - Name    =  Profil_id - Locator =  NA
## Field =  10 - Name    =  Vaissel_id - Locator =  NA

Dataset global statistics

We first establish the time amplitude of the dataset, as well as a set of colors assigned to each year.

years      = subyears = get_db_limits_year(db0)
trimesters = subtrims = seq(1,4)
colyears   = rg.colors(length(years))
cat(build_title("The dataset period is:",time2date(get_db_limits_time(db0))))
## The dataset period is: (1995-01-07 => 2016-11-29)

Let us get some statistics on the information available

db.stat.print(db0,funs=c("num","mini","maxi","mean"),
              names=c("Longitude","Latitude","Depth","Temperature","Conductivity","Salinity"))
##                 Number   Minimum   Maximum      Mean
## Longitude      5554585   -19.995    59.500    11.952
## Latitude       5554585    52.000    82.107    67.021
## Depth          5554585     0.000   100.000    49.501
## Temperature    5554574    -2.027    23.144     6.510
## Conductivity   5543746    -0.008    49.821    34.464
## Salinity       5543757    -0.000    51.512    34.525

Studying Temperature variable

From this point, most of the calculations will be performed based on the Temperature variable.

All Database

We display all samples focusing on the variable in a 2D aerial view, reporting the country borders. For comparison, we define a common color scale, established on the global minimum and maximum values (var_scale0). Note that all samples from all years and all depths are displayed (slow operation).

var = "Temperature"
colors.temp = rg.colors(rank=1)
var_scale0 = get_db_limits_var(db0,var)
display_var(db0, var = var, colors = colors.temp, title = var, filename = var)

Different Projections

We can also benefit from different projections as demonstrated next.

var = "Temperature"
filename = paste0(var," - Projection Gnomonic")
projec.define(projection="gnomonic")
## 
## Parameters for Projection
## -------------------------
## Projection is switched ON
## Use 'projec.define' to modify previous values
display_var(db0, var = var, colors = colors.temp, title = filename, filename = filename)

filename = paste0(var," - Projection Mecca[10]")
projec.define(projection="mecca",parameters=10)
## 
## Parameters for Projection
## -------------------------
## Projection is switched ON
## Use 'projec.define' to modify previous values
display_var(db0, var = var, colors = colors.temp, title = filename, filename = filename)

projec.toggle(0)

Year Campaign

We display Temperature values for the year 1995 only.

var = "Temperature"

# Comment the following line if you want to to display all years
subyears = years[1]

# Loop on the years to be displayed
for (year in subyears) 
{
  date_lim  = create_limits_date(year)
  db        = apply_sel(db0, date_lim=date_lim)
  filename  = paste0(var,"_Year_",year)
  title     = build_title(var, date_lim)
  display_var(db, var, var_scale = var_scale0, colors = colors.temp, 
              title = title, filename = filename)
}

Year/Trimester at 20m Depth

Display data for the first trimester of the year 1995 focusing on data at 20m depth exactly.

var = "Temperature"

# Select samples at depth 20m and set the color scale
db1       = apply_sel(db0, depth_lim=c(20,20), compress = TRUE)
var_scale = get_db_limits_var(db1,var)

# Comment the following lines if you want to to display all trimesters / years
subyears = years[1]
subtrims = trimesters[1]

# Loop on the years / trimesters to be displayed
for (year in subyears) 
  for (itri in subtrims)
  {
    date_lim  = create_limits_date(year, trimester=itri)
    db        = apply_sel(db1,date_lim=date_lim)
    filename  = paste0(var,"_Trim_",year,"_T",itri)
    title     = build_title(var, date_lim)
    display_var(db, var,  var_scale = var_scale, colors = colors.temp,
                title = title, filename = filename)
  }