---
title: Working with Volumes and Vectors
output: rmarkdown::html_vignette
vignette: |
  %\VignetteIndexEntry{Working with Volumes and Vectors}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
params:
  family: red
  preset: homage
css: albers.css
resource_files:
- albers.css
- albers.js
includes:
  in_header: |-
    <script src="albers.js"></script>
    <script>document.addEventListener('DOMContentLoaded',function(){document.body.classList.add('palette-red');});</script>

---

```{r setup, include = FALSE}
if (requireNamespace("ggplot2", quietly = TRUE) && requireNamespace("albersdown", quietly = TRUE)) ggplot2::theme_set(albersdown::theme_albers(family = params$family, preset = params$preset))
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  message = FALSE,
  warning = FALSE
)
suppressPackageStartupMessages(library(neuroim2))
```

`neuroim2` has two basic data containers you reach for constantly:
`NeuroVol` for one 3D image and `NeuroVec` for a 4D stack of volumes. This
article shows how to move between them, inspect geometry, and extract the
small pieces you actually analyze.

## Start with a single volume

Use `read_vol()` when the file is a 3D image or when you want one volume at a
time.

```{r read-volume}
vol_file <- system.file("extdata", "global_mask2.nii.gz", package = "neuroim2")
vol <- read_vol(vol_file)

dim(vol)
spacing(vol)
origin(vol)
```

```{r volume-check}
stopifnot(length(dim(vol)) == 3L)
stopifnot(all(spacing(vol) > 0))
```

The returned object behaves like an array, but it also carries a `NeuroSpace`
so you can keep voxel values tied to real spatial coordinates.

```{r inspect-space}
space(vol)
vol[1, 1, 1]
```

## Read a 4D time series

Use `read_vec()` when the file contains a time series or another stack of
aligned volumes.

```{r read-vector}
vec_file <- system.file("extdata", "global_mask_v4.nii", package = "neuroim2")
vec <- read_vec(vec_file)

dim(vec)
vec
```

```{r vector-check}
stopifnot(length(dim(vec)) == 4L)
stopifnot(dim(vec)[4] > 1L)
```

A `NeuroVec` still has one shared 3D spatial frame. The fourth dimension is
the axis you usually split, subset, and summarize.

## Pull out the volumes you need

`sub_vector()` is the direct way to keep only selected timepoints.

```{r subset-vector}
first_two <- sub_vector(vec, 1:2)

dim(first_two)
```

```{r subset-check}
stopifnot(dim(first_two)[4] == 2L)
```

## Convert to a matrix when you need table-like work

`as.matrix()` is useful when you want one row per voxel and one column per
timepoint.

```{r matrix-view}
mat <- as.matrix(vec)

dim(mat)
mat[1:4, 1:2]
```

```{r matrix-check}
stopifnot(nrow(mat) == prod(dim(vec)[1:3]))
stopifnot(ncol(mat) == dim(vec)[4])
```

That representation is convenient for temporal summaries, scaling, or passing
the data into general R modeling code.

## Extract one voxel or one ROI time series

For a single voxel, use `series()`.

```{r voxel-series}
voxel_ts <- series(vec, 12, 12, 12)
voxel_ts
```

For a spatial region, create an ROI and use `series_roi()`.

```{r roi-series}
roi <- spherical_roi(drop(first_two[[1]]), c(12, 12, 12), radius = 6)
roi_ts <- series_roi(vec, roi)

dim(values(roi_ts))
```

```{r roi-check}
stopifnot(length(roi) > 0L)
stopifnot(nrow(values(roi_ts)) == dim(vec)[4])
```

## Switch to sparse storage when the support matters

If only a subset of voxels should be considered present, read through a mask
and keep the result sparse.

```{r sparse-vector}
mask_vol <- read_vol(vec_file) > 0
sparse_vec <- read_vec(vec_file, mask = mask_vol)

class(sparse_vec)
dim(sparse_vec)
```

```{r sparse-check}
stopifnot(inherits(sparse_vec, "SparseNeuroVec"))
stopifnot(identical(dim(sparse_vec), dim(vec)))
```

This is the right move when absent voxels should be treated as missing support,
not as literal zeros.

## Where to go next

- `vignette("ChoosingBackends")` for storage tradeoffs
- `vignette("coordinate-systems")` for affine and orientation details
- `vignette("Resampling")` for grid changes and reorientation
- `vignette("AnalysisWorkflows")` for ROIs, searchlights, and split-map-reduce patterns
