Trying to add the ability to write.

.Rbuildignore

@@ -6,3 +6,4 @@ dim_orders.R
 test_cifti.R
 ^appveyor\.yml$
 docs/
+R/write_cifti.R

DESCRIPTION

@@ -2,7 +2,7 @@ Package: cifti
 Type: Package
 Title: Toolbox for Connectivity Informatics Technology Initiative ('CIFTI') 
     Files
-Version: 0.4.4
+Version: 0.4.5
 Authors@R: person(given = "John", family = "Muschelli", 
     role = c("aut", "cre"), email = "muschellij2@gmail.com")
 Maintainer: John Muschelli <muschellij2@gmail.com>

NAMESPACE

@@ -33,6 +33,7 @@ importFrom(gifti,readgii)
 importFrom(oro.nifti,nifti)
 importFrom(utils,download.file)
 importFrom(utils,unzip)
+importFrom(xml2,as_list)
 importFrom(xml2,read_xml)
 importFrom(xml2,xml_attrs)
 importFrom(xml2,xml_find_all)

R/get_cifti_type.R

@@ -6,6 +6,7 @@
 #' @return List of output from each type
 #' @export
 #' @importFrom R.utils gunzip
+#' @importFrom xml2 as_list
 get_cifti_type = function(
   fname,
   type = c("Volume", "Surface",
@@ -27,4 +28,25 @@ get_cifti_type = function(
   res = res[L > 0]
   return(res)
 
+}
+
+#' @rdname get_cifti_type
+cifti_as_list = function(
+  fname,
+  type = c("Volume", "Surface",
+           "Parcel",  "NamedMap",
+           "BrainModel")) {
+  type = match.arg(type, several.ok = TRUE)
+
+  nodes = matrix_ind_map_nodes(fname)
+  res = lapply(type, function(x) {
+    xml_find_all(nodes, paste0("./", x))
+  })
+  names(res) = type
+  L = sapply(res, length)
+  res = res[L > 0]
+  res = lapply(res, xml2::as_list)
+
+  return(res)
+
 }

R/nifti_2_hdr.R

@@ -7,7 +7,8 @@
 #' @return Object of class \code{nifti}
 #' @export
 #' @note The \code{unused_str} part of the header is not returned, but is an
-#' empty string of 15 characters
+#' empty string of 15 characters.  This code was adapted by
+#' the \code{oro.nifti} package
 #' @importFrom oro.nifti nifti
 nifti_2_hdr = function(fname, verbose = FALSE, warn = -1) {
   ## Open appropriate file

R/read_cifti.R

@@ -23,19 +23,22 @@ read_cifti = function(fname,
   #########################################
   # Dropping empty dimensions
   #########################################
+  orig_dim = dim(data)
   if (drop_data) {
     data = drop(data)
   }
   attr(data, "drop") = drop_data
+  attr(data, "orig_dim") = orig_dim
 
   #########################################
   # Transposing data
   #########################################
   ddata = dim(data)
-  if (length(ddata) == 2) {
-    if (trans_data) {
+  if (trans_data) {
+    if (length(ddata) == 2) {
       data = t(data)
     } else {
+      trans_data = FALSE
       warning("Dimensions of the data > 2, so no transposing done!")
     }
   }
@@ -44,6 +47,10 @@ read_cifti = function(fname,
   res$data = data
   hdr = nifti_2_hdr(fname)
   res$hdr = hdr
+  res$filename = fname
+  res$drop_data = drop_data
+  res$trans_data = trans_data
+
   class(res) = "cifti"
 
   return(res)

R/write_cifti.R

@@ -0,0 +1,166 @@
+write_nifti_2_hdr = function(hdr, filename = tempfile()) {
+
+  fid = file(filename, open = "wb")
+  hdr_size = as.integer(oro.nifti::sizeof_hdr(hdr))
+  if (hdr_size != 540L) {
+    close(fid)
+    stop("Header size is not 540 - likely not a NIfTI-2 or CIFTI file!")
+  }
+
+  validate_length = function(name, len) {
+    object = slot(hdr, name = name)
+    if (length(object) != len) {
+      stop(paste0("The ", name, " header field is not of length ", len))
+    }
+  }
+
+  writeBin(as.integer(hdr@"sizeof_hdr"), fid, size=4)
+  writeChar(hdr@magic, fid, nchars = 8, eos = NULL)
+  writeBin(as.integer(hdr@"datatype"), fid, size=2)
+  writeBin(as.integer(hdr@"bitpix"), fid, size=2)
+
+  # worked for int64
+  validate_length("dim_", 8)
+  writeBin(as.integer(hdr@"dim_"), fid, size = 8)
+  writeBin(as.double(hdr@"intent_p1"), fid, size=8)
+  writeBin(as.double(hdr@"intent_p2"), fid, size=8)
+  writeBin(as.double(hdr@"intent_p3"), fid, size=8)
+  writeBin(as.double(hdr@"pixdim"), fid, size=8)
+
+  if (hdr@"vox_offset" < 544) {
+    warning("vox_offset seems off!")
+  }
+  writeBin(as.integer(hdr@"vox_offset"), fid, size=8)
+
+  writeBin(as.double(hdr@"scl_slope"), fid, size=8)
+  writeBin(as.double(hdr@"scl_inter"), fid, size=8)
+  writeBin(as.double(hdr@"cal_max"), fid, size=8)
+  writeBin(as.double(hdr@"cal_min"), fid, size=8)
+
+
+  writeBin(as.double(hdr@"slice_duration"), fid, size=8)
+  writeBin(as.double(hdr@"toffset"), fid, size=8)
+
+
+  writeBin(as.integer(hdr@"slice_start"), fid, size=8)
+  writeBin(as.integer(hdr@"slice_end"), fid, size=8)
+
+  writeChar(hdr@descrip, fid, nchars = 80, eos = NULL)
+  writeChar(hdr@aux_file, fid, nchars = 24, eos = NULL)
+
+  writeBin(as.integer(hdr@"qform_code"), fid, size=4)
+  writeBin(as.integer(hdr@"sform_code"), fid, size=4)
+
+  writeBin(as.double(hdr@"quatern_b"), fid, size=8)
+  writeBin(as.double(hdr@"quatern_c"), fid, size=8)
+  writeBin(as.double(hdr@"quatern_d"), fid, size=8)
+  writeBin(as.double(hdr@"qoffset_x"), fid, size=8)
+  writeBin(as.double(hdr@"qoffset_y"), fid, size=8)
+  writeBin(as.double(hdr@"qoffset_z"), fid, size=8)
+
+  validate_length("srow_x", 4)
+  validate_length("srow_y", 4)
+  validate_length("srow_z", 4)
+
+  writeBin(as.double(hdr@"srow_x"), fid, size=8)
+  writeBin(as.double(hdr@"srow_y"), fid, size=8)
+  writeBin(as.double(hdr@"srow_z"), fid, size=8)
+
+  writeBin(as.integer(hdr@"slice_code"), fid, size=4)
+  writeBin(as.integer(hdr@"xyzt_units"), fid, size=4)
+  writeBin(as.integer(hdr@"intent_code"), fid, size=4)
+
+  writeChar(hdr@intent_name, fid, nchars = 16, eos = NULL)
+  writeChar(hdr@dim_info, fid, nchars = 1, eos = NULL)
+
+  suppressWarnings({
+    writeChar("", fid, nchars = 15, eos = NULL)
+  })
+
+  nhdr = seek(fid)
+  stopifnot(nhdr == hdr_size)
+  return(list(fid = fid, filename = filename))
+
+}
+
+
+
+
+
+write_cifti = function(res) {
+
+  filename = tempfile()
+  hdr = res$hdr
+  L = write_nifti_2_hdr(hdr, filename = filename)
+  fid = L$fid
+
+  data = res$data
+
+  # reversing things that happened with read_cifti
+  ddata = dim(data)
+  trans_data = attr(data, "trans")
+  if (is.null(trans_data)) {
+    trans_data = FALSE
+  }
+  if (trans_data) {
+    if (length(ddata) == 2) {
+      data = t(data)
+    } else {
+      trans_data = FALSE
+      warning("Dimensions of the data > 2, so no transposing done!")
+    }
+  }
+
+  # reversing things that happened with read_cifti
+  drop_data = attr(data, "drop")
+  if (is.null(drop_data)) {
+    drop_data = FALSE
+  }
+  orig_dim = attr(data, "orig_dim")
+  if (drop_data) {
+    data = array(data, dim = orig_dim)
+  }
+
+
+  seek(fid, where = hdr@vox_offset, origin = "start");
+
+  dtype = as.character(hdr@datatype)
+  what_func = switch(
+    dtype,
+    "2" = as.integer,
+    "4" = as.numeric,
+    "8" = as.integer,
+    "16" = as.numeric,
+    "64" = as.double,
+    "512" = as.numeric,
+    "768" = as.integer
+  )
+
+  size = switch(
+    dtype,
+    "2" = 1,
+    "4" = 2,
+    "8" = 4,
+    "16" = 4,
+    "64" = 8,
+    "512" = 2,
+    "768" = 4
+  )
+  if (is.null(what_func) || is.null(size)) {
+    stop("Unsupported data type indicated by NIfTI-2 header!")
+  }
+
+  vals = c(data)
+  img_dim = hdr@dim_
+  n_items = prod(img_dim[2:length(img_dim)])
+  if (n_items > length(vals)) {
+    stop("Not all CIFTI data read!")
+  }
+  if (n_items < length(vals)) {
+    stop("More data read than header indicated - header or data problem!")
+  }
+
+  cifti_dim = img_dim[6:8]
+  vals = array(vals, dim = cifti_dim)
+
+}