Nothing
##
##
## Copyright (c) 2009-2014 Brandon Whitcher and Volker Schmid
## All rights reserved.
##
## Redistribution and use in source and binary forms, with or without
## modification, are permitted provided that the following conditions are
## met:
##
## * Redistributions of source code must retain the above copyright
## notice, this list of conditions and the following disclaimer.
## * Redistributions in binary form must reproduce the above
## copyright notice, this list of conditions and the following
## disclaimer in the documentation and/or other materials provided
## with the distribution.
## * The names of the authors may not be used to endorse or promote
## products derived from this software without specific prior
## written permission.
##
## THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
## "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
## LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
## A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
## HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
## SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
## LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
## DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
## THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
## (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
## OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
##
## $Id: readS4.R 332 2010-01-29 16:54:07Z bjw34032 $
##
## Sub-routines
.readCharWithEmbeddedNuls <- function(fid, n, to="UTF-8") {
txt <- readBin(fid, "raw", n)
iconv(rawToChar(txt[txt != as.raw(0)]), to=to)
}
##
## readNIfTI() is a convient interface for the user
##
readNIfTI <- function(fname, verbose=FALSE, warn=-1, reorient=TRUE,
call=NULL) {
if (is.null(call)) {
call <- match.call()
}
## Warnings?
oldwarn <- getOption("warn")
options(warn=warn)
if (verbose) {
cat(paste(" fname =", fname), fill=TRUE)
}
## Separate path and file name for more robust extension stripping
pathVector <- unlist(strsplit(fname, "/"))
file.name <- pathVector[length(pathVector)]
path <- paste(pathVector[-length(pathVector)], collapse="/")
if (length(pathVector) > 1) {
fname <- paste(path, file.name, sep="/")
} else {
fname <- file.name
}
## Strip any extensions
fname <- sub("\\.gz$", "", fname)
fname <- sub("\\.nii$", "", fname)
fname <- sub("\\.hdr$", "", fname)
fname <- sub("\\.img$", "", fname)
## Add all possible file extensions
nii <- paste(fname, "nii", sep=".")
niigz <- paste(fname, "nii.gz", sep=".")
hdr <- paste(fname, "hdr", sep=".")
hdrgz <- paste(fname, "hdr.gz", sep=".")
img <- paste(fname, "img", sep=".")
imggz <- paste(fname, "img.gz", sep=".")
## Check all possible file extensions
if (file.exists(niigz)) {
## If compressed file exists, then upload!
if (verbose) {
cat(paste(" files =", niigz), fill=TRUE)
}
nim <- .read.nifti.content(fname, onefile=TRUE, gzipped=TRUE,
verbose=verbose, warn=warn, reorient=reorient,
call=call)
} else {
if (file.exists(nii)) {
## If uncompressed file exists, then upload!
if (verbose) {
cat(paste(" files =", nii), fill=TRUE)
}
nim <- .read.nifti.content(fname, onefile=TRUE, gzipped=FALSE,
verbose=verbose, warn=warn, reorient=reorient,
call=call)
} else {
if (file.exists(hdrgz) && file.exists(imggz)) {
## If compressed files exist, then upload!
if (verbose) {
cat(paste(" files =", hdrgz, "and", imggz), fill=TRUE)
}
nim <- .read.nifti.content(fname, onefile=FALSE, gzipped=TRUE,
verbose=verbose, warn=warn,
reorient=reorient, call=call)
} else {
## If uncompressed files exist, then upload!
if (file.exists(hdr) && file.exists(img)) {
if (verbose) {
cat(paste(" files =", hdr, "and", img), fill=TRUE)
}
nim <- .read.nifti.content(fname, onefile=FALSE, gzipped=FALSE,
verbose=verbose, warn=warn,
reorient=reorient, call=call)
} else {
stop("File(s) not found!")
}
}
}
}
options(warn=oldwarn)
return(nim)
}
############################################################################
############################################################################
############################################################################
.read.nifti.content <- function(fname, onefile=TRUE, gzipped=TRUE,
verbose=FALSE, warn=-1, reorient=FALSE,
call=NULL) {
## Open appropriate file
if (gzipped) {
suffix <- ifelse(onefile, "nii.gz", "hdr.gz")
fname <- paste(fname, suffix, sep=".")
fid <- gzfile(fname, "rb")
if (verbose) {
cat(" nii =", fname, fill=TRUE)
}
} else {
suffix <- ifelse(onefile, "nii", "hdr")
fname <- paste(fname, suffix, sep=".")
fid <- file(fname, "rb")
if (verbose) {
cat(" hdr =", fname, fill=TRUE)
}
}
## Warnings?
oldwarn <- getOption("warn")
options(warn=warn)
## Test for endian properties
endian <- .Platform$endian
sizeof.hdr <- readBin(fid, integer(), size=4, endian=endian)
if (sizeof.hdr != 348) {
close(fid)
endian <- "swap"
if (gzipped) {
fid <- gzfile(fname, "rb")
} else {
fid <- file(fname, "rb")
}
sizeof.hdr <- readBin(fid, integer(), size=4, endian=endian)
if (verbose) {
cat(" ENDIAN = swap", fill=TRUE)
}
}
## Construct S4 object
nim <- nifti()
nim@"sizeof_hdr" <- sizeof.hdr
nim@"data_type" <- .readCharWithEmbeddedNuls(fid, n=10)
nim@"db_name" <- .readCharWithEmbeddedNuls(fid, n=18)
nim@"extents" <- readBin(fid, integer(), size=4, endian=endian)
nim@"session_error" <- readBin(fid, integer(), size=2, endian=endian)
nim@"regular" <- .readCharWithEmbeddedNuls(fid, n=1)
nim@"dim_info" <- .readCharWithEmbeddedNuls(fid, n=1)
nim@"dim_" <- readBin(fid, integer(), 8, size=2, endian=endian)
nim@"intent_p1" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"intent_p2" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"intent_p3" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"intent_code" <- readBin(fid, integer(), size=2, endian=endian)
nim@"datatype" <- readBin(fid, integer(), size=2, endian=endian)
nim@"bitpix" <- readBin(fid, integer(), size=2, endian=endian)
nim@"slice_start" <- readBin(fid, integer(), size=2, endian=endian)
nim@"pixdim" <- readBin(fid, numeric(), 8, size=4, endian=endian)
nim@"vox_offset" <- readBin(fid, numeric(), size=4, endian=endian)
if (verbose) {
cat(" vox_offset =", nim@"vox_offset", fill=TRUE)
}
nim@"scl_slope" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"scl_inter" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"slice_end" <- readBin(fid, integer(), size=2, endian=endian)
nim@"slice_code" <- readBin(fid, integer(), size=1, signed=FALSE,
endian=endian)
nim@"xyzt_units" <- readBin(fid, integer(), size=1, signed=FALSE,
endian=endian)
nim@"cal_max" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"cal_min" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"slice_duration" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"toffset" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"glmax" <- readBin(fid, integer(), size=4, endian=endian)
nim@"glmin" <- readBin(fid, integer(), size=4, endian=endian)
nim@"descrip" <- .readCharWithEmbeddedNuls(fid, n=80)
nim@"aux_file" <- .readCharWithEmbeddedNuls(fid, n=24)
nim@"qform_code" <- readBin(fid, integer(), size=2, endian=endian)
nim@"sform_code" <- readBin(fid, integer(), size=2, endian=endian)
nim@"quatern_b" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"quatern_c" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"quatern_d" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"qoffset_x" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"qoffset_y" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"qoffset_z" <- readBin(fid, numeric(), size=4, endian=endian)
nim@"srow_x" <- readBin(fid, numeric(), 4, size=4, endian=endian)
nim@"srow_y" <- readBin(fid, numeric(), 4, size=4, endian=endian)
nim@"srow_z" <- readBin(fid, numeric(), 4, size=4, endian=endian)
nim@"intent_name" <- .readCharWithEmbeddedNuls(fid, n=16)
nim@"magic" <- .readCharWithEmbeddedNuls(fid, n=4)
## To flag such a struct as being conformant to the NIFTI-1 spec,
## the last 4 bytes of the header must be either the C String "ni1"
## or "n+1"; in hexadecimal, the 4 bytes [6E 69 31 00] or [6E 2B 31
## 00] (in any future version of this format, the 1 will be upgraded
## to 2, etc.). Normally, such a "magic number" or flag goes at the
## start of the file, but trying to avoid clobbering widely-used
## ANALYZE 7.5 fields led to putting this marker last. However,
## recall that "the last shall be first" (Matthew 20:16).
if (onefile) {
if (nim@"magic" != "n+1") {
stop("This is not a one-file NIfTI format")
}
nim@"extender" <- readBin(fid, integer(), 4, size=1, signed=FALSE,
endian=endian)
## If extension[0] is nonzero, it indicates that extended header
## information is present in the bytes following the extension
## array. In a .nii file, this extended header data is before the
## image data (and vox_offset must be set correctly to allow for
## this). In a .hdr file, this extended data follows extension and
## proceeds (potentially) to the end of the file.
##
if (nim@"extender"[1] > 0 || nim@"vox_offset" > 352) {
if (verbose) {
cat(" niftiExtension detected!", fill=TRUE)
}
if (!is(nim, "niftiExtension")) {
nim <- as(nim, "niftiExtension")
}
while (seek(fid) < nim@"vox_offset") {
if (verbose) {
cat(" seek(fid) =", seek(fid), fill=TRUE)
}
nimextsec <- new("niftiExtensionSection")
nimextsec@esize <- readBin(fid, integer(), size=4, endian=endian)
nimextsec@ecode <- readBin(fid, integer(), size=4, endian=endian)
nimextsec@edata <- .readCharWithEmbeddedNuls(fid, n=nimextsec@esize-8)
nim@extensions <- append(nim@extensions, nimextsec)
}
if (seek(fid) > nim@"vox_offset") {
stop("-- extension size (esize) has overshot voxel offset --")
}
}
}
if (verbose) {
cat(" seek(fid) =", seek(fid), fill=TRUE)
}
dims <- 2:(1+nim@"dim_"[1])
n <- prod(nim@"dim_"[dims])
if (! onefile) {
if (nim@"magic" != "ni1") {
stop("This is not a two-file NIfTI format")
}
close(fid)
fname <- sub("\\.hdr$", "\\.img", fname)
if (gzipped) {
fid <- gzfile(fname, "rb")
} else {
fid <- file(fname, "rb")
}
## seek(fid, nim@"vox_offset") # not necessary for two-file format
}
data <-
switch(as.character(nim@"datatype"),
"2" = readBin(fid, integer(), n, nim@"bitpix"/8, signed=FALSE,
endian=endian),
"4" = readBin(fid, integer(), n, nim@"bitpix"/8, endian=endian),
"8" = readBin(fid, integer(), n, nim@"bitpix"/8, endian=endian),
"16" = readBin(fid, double(), n, nim@"bitpix"/8, endian=endian),
"64" = readBin(fid, double(), n, nim@"bitpix"/8, endian=endian),
"512" = readBin(fid, integer(), n, nim@"bitpix"/8, signed=FALSE,
endian=endian),
"768" = readBin(fid, integer(), n, nim@"bitpix"/8, signed=FALSE,
endian=endian),
stop(paste("Data type", nim@"datatype", "unsupported in", fname))
)
close(fid)
## WARNING to the user
if (nim@"scl_slope" != 0) {
warning(paste("scl_slope =", nim@"scl_slope", "and data must be rescaled."))
data <- data * nim@"scl_slope" + nim@"scl_inter"
}
##
## THE SLOW BIT FOLLOWS
##
## 3D IMAGE (VOLUME) ORIENTATION AND LOCATION IN SPACE:
## There are 3 different methods by which continuous coordinates can
## attached to voxels. The discussion below emphasizes 3D volumes,
## and the continuous coordinates are referred to as (x,y,z). The
## voxel index coordinates (i.e., the array indexes) are referred to
## as (i,j,k), with valid ranges:
## i = 0 .. dim[1]-1
## j = 0 .. dim[2]-1 (if dim[0] >= 2)
## k = 0 .. dim[3]-1 (if dim[0] >= 3)
## The (x,y,z) coordinates refer to the CENTER of a voxel. In
## methods 2 and 3, the (x,y,z) axes refer to a subject-based
## coordinate system, with
## +x = Right +y = Anterior +z = Superior.
## This is a right-handed coordinate system. However, the exact
## direction these axes point with respect to the subject depends on
## qform_code (Method 2) and sform_code (Method 3).
##
if (reorient) {
nim@.Data <- reorient(nim, data, verbose=verbose)
nim@"reoriented" <- TRUE
} else {
nim@.Data <- array(data, nim@"dim_"[dims])
}
## Warnings?
options(warn=oldwarn)
## Check validity
validNIfTI <- getValidity(getClassDef("nifti"))
validNIfTI(nim)
if (getOption("niftiAuditTrail")) {
if (is.null(call)) {
call <- match.call()
}
nim <- niftiExtensionToAuditTrail(nim, workingDirectory=getwd(),
filename=fname, call=call)
}
return(nim)
}
############################################################################
## readANALYZE() is a convenient interface for the user
############################################################################
readANALYZE <- function(fname, SPM=FALSE, verbose=FALSE, warn=-1) {
## Warnings?
oldwarn <- getOption("warn")
options(warn=warn)
if (verbose) {
cat(paste(" fname =", fname), fill=TRUE)
}
## Separate path and file name for more robust extension stripping
pathVector <- unlist(strsplit(fname, "/"))
file.name <- pathVector[length(pathVector)]
path <- paste(pathVector[-length(pathVector)], collapse="/")
if (length(pathVector) > 1) {
fname <- paste(path, file.name, sep="/")
} else {
fname <- file.name
}
## Strip any extensions
fname <- sub("\\.gz$", "", fname)
fname <- sub("\\.hdr$", "", fname)
fname <- sub("\\.img$", "", fname)
if (! (file.exists(paste(fname, "hdr", sep=".")) &&
file.exists(paste(fname, "img", sep="."))) &&
! (file.exists(paste(fname, "hdr.gz", sep=".")) &&
file.exists(paste(fname, "img.gz", sep=".")))) {
stop("File(s) not found!")
}
## If uncompressed files exist, then upload!
if (file.exists(paste(fname, "hdr", sep=".")) &&
file.exists(paste(fname, "img", sep="."))) {
if (verbose) {
cat(paste(" files = ", fname, ".{hdr,img}", sep=""), fill=TRUE)
}
aim <- .read.analyze.content(fname, gzipped=FALSE, SPM=SPM,
verbose=verbose, warn=warn)
options(warn=oldwarn)
return(aim)
}
## If compressed files exist, then upload!
if (file.exists(paste(fname, "hdr.gz", sep=".")) &&
file.exists(paste(fname, "img.gz", sep="."))) {
if (verbose) {
cat(paste(" files = ", fname, ".{hdr.gz,img.gz}", sep=""), fill=TRUE)
}
aim <- .read.analyze.content(fname, gzipped=TRUE, SPM=SPM,
verbose=verbose, warn=warn)
options(warn=oldwarn)
return(aim)
}
invisible()
}
############################################################################
############################################################################
############################################################################
.read.analyze.content <- function(fname, gzipped=TRUE, SPM=FALSE,
verbose=FALSE, warn=-1) {
## Open header file
if (gzipped) {
fname <- paste(fname, "hdr.gz", sep=".")
fid <- gzfile(fname, "rb")
} else {
fname <- paste(fname, "hdr", sep=".")
fid <- file(fname, "rb")
}
if (verbose) {
cat(" hdr =", fname, fill=TRUE)
}
## Warnings?
oldwarn <- getOption("warn")
options(warn=warn)
## Test for endian properties
endian <- .Platform$endian
sizeof.hdr <- readBin(fid, integer(), size=4, endian=endian)
if (sizeof.hdr != 348) {
close(fid)
endian <- "swap"
if (gzipped) {
fid <- gzfile(fname, "rb")
} else {
fid <- file(fname, "rb")
}
sizeof.hdr <- readBin(fid, integer(), size=4, endian=endian)
}
## Construct S4 object
aim <- new("anlz")
aim@"sizeof_hdr" <- sizeof.hdr
aim@"data_type" <- .readCharWithEmbeddedNuls(fid, n=10)
aim@"db_name" <- .readCharWithEmbeddedNuls(fid, n=18)
aim@"extents" <- readBin(fid, integer(), size=4, endian=endian)
aim@"session_error" <- readBin(fid, integer(), size=2, endian=endian)
aim@"regular" <- .readCharWithEmbeddedNuls(fid, n=1)
aim@"hkey_un0" <- .readCharWithEmbeddedNuls(fid, n=1)
aim@"dim_" <- readBin(fid, integer(), 8, size=2, endian=endian)
aim@"vox_units" <- .readCharWithEmbeddedNuls(fid, n=4)
aim@"cal_units" <- .readCharWithEmbeddedNuls(fid, n=8)
aim@"unused1" <- readBin(fid, integer(), size=2, endian=endian)
aim@"datatype" <- readBin(fid, integer(), size=2, endian=endian)
aim@"bitpix" <- readBin(fid, integer(), size=2, endian=endian)
aim@"dim_un0" <- readBin(fid, integer(), size=2, endian=endian)
aim@"pixdim" <- readBin(fid, numeric(), 8, size=4, endian=endian)
aim@"vox_offset" <- readBin(fid, numeric(), size=4, endian=endian)
aim@"funused1" <- readBin(fid, numeric(), size=4, endian=endian)
## SPM has used the ANALYZE 7.5 funused1 field as a scaling factor
if (aim@"funused1" != 0) {
warning(paste("funused1 =", aim@"funused1", "and data must be rescaled."))
}
##
aim@"funused2" <- readBin(fid, numeric(), size=4, endian=endian)
## Maybe I'm paranoid, but let's check funused2 in case it is an intercept
if (aim@"funused2" != 0) {
warning(paste("funused2 =", aim@"funused2", "and data must be translated."))
}
##
aim@"funused3" <- readBin(fid, numeric(), size=4, endian=endian)
aim@"cal_max" <- readBin(fid, numeric(), size=4, endian=endian)
aim@"cal_min" <- readBin(fid, numeric(), size=4, endian=endian)
aim@"compressed" <- readBin(fid, numeric(), size=4, endian=endian)
aim@"verified" <- readBin(fid, numeric(), size=4, endian=endian)
aim@"glmax" <- readBin(fid, integer(), size=4, endian=endian)
aim@"glmin" <- readBin(fid, integer(), size=4, endian=endian)
aim@"descrip" <- .readCharWithEmbeddedNuls(fid, n=80)
aim@"aux_file" <- .readCharWithEmbeddedNuls(fid, n=24)
aim@"orient" <- .readCharWithEmbeddedNuls(fid, n=1)
aim@"origin" <- readBin(fid, integer(), 5, size=2, endian=endian) # .readCharWithEmbeddedNuls(fid, n=10)
aim@"generated" <- .readCharWithEmbeddedNuls(fid, n=10)
aim@"scannum" <- .readCharWithEmbeddedNuls(fid, n=10)
aim@"patient_id" <- .readCharWithEmbeddedNuls(fid, n=10)
aim@"exp_date" <- .readCharWithEmbeddedNuls(fid, n=10)
aim@"exp_time" <- .readCharWithEmbeddedNuls(fid, n=10)
aim@"hist_un0" <- .readCharWithEmbeddedNuls(fid, n=3)
aim@"views" <- readBin(fid, integer(), size=4, endian=endian)
aim@"vols_added" <- readBin(fid, integer(), size=4, endian=endian)
aim@"start_field" <- readBin(fid, integer(), size=4, endian=endian)
aim@"field_skip" <- readBin(fid, integer(), size=4, endian=endian)
aim@"omax" <- readBin(fid, integer(), size=4, endian=endian)
aim@"omin" <- readBin(fid, integer(), size=4, endian=endian)
aim@"smax" <- readBin(fid, integer(), size=4, endian=endian)
magic <- readBin(fid, raw(), n=4, endian=endian)
aim@"smin" <- readBin(magic, integer(), size=4, endian=endian)
## Test for "magic" field (should not exist)
if (rawToChar(magic) == "ni1") { # now its actually NIfTI two-file format
stop("This is in two-file NIfTI format, please use readNIfTI")
}
close(fid)
## Open image file
if (gzipped) {
fname <- sub("\\.hdr\\.gz$", "\\.img\\.gz", fname) # paste(fname, ".img.gz", sep=".")
fid <- gzfile(fname, "rb")
} else {
fname <- sub("\\.hdr$", "\\.img", fname) # paste(fname, "img", sep=".")
fid <- file(fname, "rb")
}
if (verbose) {
cat(" img =", fname, fill=TRUE)
}
n <- prod(aim@"dim_"[2:5])
data <- switch(as.character(aim@"datatype"),
"1" = readBin(fid, integer(), n, aim@"bitpix"/8,
signed=FALSE, endian=endian),
"2" = readBin(fid, integer(), n, aim@"bitpix"/8,
signed=FALSE, endian=endian),
"4" = readBin(fid, integer(), n, aim@"bitpix"/8,
endian=endian),
"8" = readBin(fid, integer(), n, aim@"bitpix"/8,
endian=endian),
"16" = readBin(fid, numeric(), n, aim@"bitpix"/8,
endian=endian),
"64" = readBin(fid, double(), n, aim@"bitpix"/8,
endian=endian),
stop(paste("Data type ", aim@"datatype", " (",
convert.datatype.anlz(aim@"datatype"),
") unsupported in", fname, sep="")))
close(fid)
dims <- 2:(1+aim@"dim_"[1])
if (SPM) {
if (verbose) {
cat(" SPM format has been specified and data re-scaled.", fill=TRUE)
}
aim@.Data <- array(aim@"funused1" * data, aim@"dim_"[dims])
} else {
aim@.Data <- array(data, aim@"dim_"[dims])
}
## Warnings?
options(warn=oldwarn)
## Check validity
validANALYZE <- getValidity(getClassDef("anlz"))
validANALYZE(aim)
return(aim)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.