inst/third_party/AFNIio.R
In beauchamplab/rave: R Analysis and Visualization of 'ECOG/iEEG' Data
#------------------------------------------------------------------
# Global Variables
#------------------------------------------------------------------
BATCH_MODE <<- 0 #Initialize batch mode flag to 0
R_io <<- -1
SHOW_TRC <<- FALSE
#------------------------------------------------------------------
# Functions for library loading
#------------------------------------------------------------------
find.in.path <- function(file) { #Pretty much same as first.in.path
sys_path = strsplit(Sys.getenv('PATH'),':')[[1]]
ff <- paste(sys_path,'/', file, sep='')
ff = ff[file.exists(ff)]
ff = base::normalizePath(ff)
aa <- gsub('//','/',ff[1], fixed=TRUE)
if (is.na(aa)) aa <- NULL
return(aa)
}
#Return cat's file option to format the help output
#for TXT or sphinx purposes
help.cat.file.AFNI <- function (pname=NULL, targ='TXT') {
if (is.null(pname)) {
err.AFNI("NULL name for help function")
return("")
}
if (targ == 'TXT') {
dopt = '-hdoc_2_txt'
} else if (targ == 'SPX') {
dopt = '-hdoc_2_spx'
} else if (targ == 'ASPX') {
dopt = '-hdoc_2_aspx'
} else if (targ == 'RAW') {
return("")
} else {
warn.AFNI(paste("targ", targ,"unknown. Assuming 'TXT'"))
dopt = '-hdoc_2_txt'
}
return(paste("|apsearch ", dopt ,pname,"-"))
}
#print warnings a la AFNI
prompt.AFNI <- function (str='I await', choices=c('y','n'), vals=NULL) {
if (!is.null(vals) && length(vals) != length(choices)) {
err.AFNI(paste("Have ", length(choices), "options, but",
length(vals), "values to return"))
return(0)
}
choices[1]<-toupper(choices[1])
kk<-vector(length=0)
spr <- paste(str," [",paste(choices, collapse='|'),"]:", sep='')
if (BATCH_MODE) { #Only two choices available for this one!
if (length(choices)!=2) {
err.AFNI("This one can't run in batch mode for more than 2 choices")
return(0)
}
ss <- sys.AFNI(paste(
"prompt_user -pause '", spr,"'", sep='' ))
if (ss$out == "0") return(2) #The not default
else return(1) #The default
} else {
while (length(kk) == 0) {
cat(spr)
bb <- readLines(n=1)
if (bb == '') {
kk <- 1
} else {
kk <-which(tolower(choices) == tolower(bb))
}
}
if (!is.null(vals)) {
return(vals[kk])
} else {
return(kk)
}
}
return(0)
}
set.AFNI.msg.trace <- function (vv=FALSE) { SHOW_TRC <<- vv }
warn.AFNI <- function (str='Consider yourself warned',
callstr=NULL,
newline=TRUE) {
if (is.null(callstr)) {
if (SHOW_TRC) callstr <- who.called.me(TRUE)
else callstr <- ''
}
nnn<-''
if (newline) nnn <- '\n'
if (BATCH_MODE) ff <- stderr()
else ff <- ''
cat( '\n', 'oo Warning: ', callstr,'\n ',
paste(str, collapse=''), nnn,
sep='', file = ff)
}
err.AFNI <- function (str='Danger Danger Will Robinson',
callstr=NULL,
newline=TRUE) {
if (is.null(callstr)) {
if (SHOW_TRC) callstr <- who.called.me(TRUE)
else callstr <- ''
}
nnn<-''
if (newline) nnn <- '\n'
if (BATCH_MODE) ff <- stderr()
else ff <- ''
cat( '\n', '** Error: ', callstr,'\n ',
paste(str, collapse=''), nnn,
sep='', file = ff)
}
errex.AFNI <- function (str='Alas this must end',
callstr=NULL, newline=TRUE) {
if (is.null(callstr)) {
if (SHOW_TRC) callstr <- who.called.me(TRUE)
else callstr <- ''
}
err.AFNI(str,callstr, newline)
exit.AFNI(str='\n Execution halted',stat=1)
}
note.AFNI <- function (str='May I speak frankly?',
callstr=NULL, newline=TRUE, tic=1,
trimtrace=30) {
if (is.null(callstr)) {
if (SHOW_TRC) callstr <- who.called.me(TRUE, trim=trimtrace)
else callstr <- ''
}
nnn<-''
if (newline) nnn <- '\n'
if (BATCH_MODE) ff <- stderr()
else ff <- ''
if (tic == 1) {
tm <- format(Sys.time(), " @ %H:%M:%S")
} else if (tic==2) {
tm <- format(Sys.time(), " @ %a %b %d %H:%M:%S %Y")
} else tm <- ''
cat( '\n', '++ Note: ', callstr,
sprintf('%s\n ', tm),
paste(str, collapse=''),nnn,
sep='', file = ff)
}
exit.AFNI <- function(str='The piano has been drinking.', stat=0) {
if (BATCH_MODE) {
quit(save='no', status = stat)
} else {
#note.AFNI(str)
stop(str)
}
}
#Locate and load R_io.so
set_R_io <- function() {
rio <- 0
ll <- find.in.path('R_io.so')
if (!is.null(ll)) {
dd <- try(dyn.load(ll), silent=TRUE)
if (dd[[1]]!="R_io") {
warn.AFNI(paste("Failed to load R_io.so with this error message:\n"))
dyn.load(ll)
} else {
rio <- 1
}
}
return(rio)
}
if (R_io == -1) {
R_io <<- set_R_io()
}
have_R_io <- function() {
if (R_io == 1) return(TRUE) else return(FALSE)
}
libLoad <- function(myLib) {
sucLoad <- FALSE
sucCheck <- FALSE
try(sucLoad <- library(myLib, character.only = TRUE, logical.return = TRUE))
if (sucLoad) {
print(sprintf("Package %s successfully loaded!", myLib)); sucCheck <- TRUE
} else {
if (BATCH_MODE == 1) {
err.AFNI(paste( "Need to install package ",myLib,
"\n Start an interactive R session then run:\n",
"\n update.packages(checkBuilt=TRUE, ask=FALSE) ",
"#Not mandatory, but recommended",
"\n install.packages('",myLib,"')\n",
"\n Then quit R and rerun your command.\n",
sep=''))
} else {
try(install.packages(myLib))
try(sucLoad <- library(myLib, character.only = TRUE,
logical.return = TRUE))
if (sucLoad) print(sprintf("Package %s successfully loaded...", myLib))
}
}
return(sucLoad)
}
pkgLoad <- function(pkgs) {
# install packages not already loaded:
pkgs_miss <- pkgs[which(!pkgs %in% installed.packages()[, 1])]
if (length(pkgs_miss) > 0) {
install.packages(pkgs_miss, dep=TRUE, repos='http://watson.nci.nih.gov/cran_mirror/')
}
# load packages not already loaded:
attached <- search()
attached_pkgs <- attached[grepl("package", attached)]
need_to_attach <- pkgs[which(!pkgs %in% gsub("package:", "", attached_pkgs))]
if (length(need_to_attach) > 0) {
for (i in 1:length(need_to_attach)) if(require(need_to_attach[i], character.only = TRUE))
cat('Package ', need_to_attach[i], ' loaded successfully!\n\n', sep='')
}
if (length(need_to_attach) == 0) {
message("\n ...All required packages were already loaded!\n")
}
}
#------------------------------------------------------------------
# Functions to deal with AFNI file names
#------------------------------------------------------------------
strip.extension <- function (filename, extvec=NULL, verb=0) {
n <- list()
if (is.null(extvec)) {
ff <- strsplit(filename, '\\.')[[1]]
if (length(ff) > 1) {
n$ext <- paste('.',ff[length(ff)], sep='')
n$name_noext <- paste(ff[1:length(ff)-1],collapse='.')
} else {
n$ext <- ''
n$name_noext <- filename
}
} else {
n$ext <- ''
n$name_noext <- filename
for (ex in extvec) {
patt <- paste('\\',ex,'$',collapse='', sep='')
if (length(grep(patt, filename))) {
n$ext <- ex
n$name_noext <- sub(patt,'',filename)
return(n)
}
}
}
return(n)
}
parse.name <- function (filename, extvec=NULL, verb=0) {
n <- list()
ff <- strsplit(filename, .Platform$file.sep)[[1]]
n$name <- ff[length(ff)]
n$path <- '.'
if (length(ff)>1)
n$path <- paste(ff[1:length(ff)-1],collapse=.Platform$file.sep)
n$path <- paste(n$path, .Platform$file.sep, sep="")
n2 <- strip.extension(n$name, extvec, verb)
n$ext <- n2$ext
n$name_noext <- n2$name_noext
return(n)
}
is.AFNI.1D.string <- function(t) {
if (length(grep('^1D:',t))) return(TRUE)
return(FALSE)
}
eval.AFNI.string.help <- function() {
return("
Data Strings:
-------------
You can specify input matrices and vectors in a variety of
ways. The simplest is by specifying a .1D file with all
the trimmings of column and row selectors. You can also
specify a string that gets evaluated on the fly.
For example: '1D: 1 4 8' evaluates to a vector of values 1 4 and 8.
Also, you can use R expressions such as: 'R: seq(0,10,3)'
")
}
eval.AFNI.1D.string <- function (t, verb=0, nmax=0) {
#remove 1D:
t <- sub("^1D:","",t)
#any transpose?
doTr = FALSE
if (length(grep("'$",t))) {
t<-sub("'$",'',t)
doTr = TRUE
}
#replace commas with space
t<-gsub(",",' ',t)
#remove multiple blanks
t <- deblank.string(t, middle=TRUE)
vvf <- vector(length = 0, mode="numeric")
#replace .. with : and split at 'space'
s <- strsplit(sub("..",":",t, fixed=TRUE), " ")[[1]]
#replace $ with nmax if possible
s <- gsub("[$]",as.character(nmax),s)
#Now loop and form vector of components
for (si in s) {
if (verb) cat ("working ", si, "\n")
if (length(grep(":",si))) {
sss <- strsplit(si,'[(*)]')[[1]]
if (length(sss)>1) {
se <- sss[2]
si <- strsplit(sss[1],":")[[1]]
vv <- eval(parse(text=sprintf('seq(from=%s,to=%s,by=%s)',
si[1], si[2], se) ))
} else {
vv <- eval(parse(text=si))
}
} else if (length(grep("@",si))) {
ssi = as.numeric(strsplit(si,"@")[[1]])
#cat ("ssi = ", ssi, "\n")
vv = rep(ssi[2], ssi[1])
} else {
vv = as.numeric(si)
}
#cat(si," = ",vv, "\n")
vvf <- c(vvf, vv)
#cat("vvnow = ",vvf, "\n")
}
if (doTr) return(t(vvf))
else return(vvf)
}
is.AFNI.R.string <- function(t) {
if (length(grep('^R:',t))) return(TRUE)
return(FALSE)
}
eval.AFNI.R.string <- function (t) {
t <- sub("^R:","",t)
return(eval(parse(text=t)))
}
eval.AFNI.string <- function (t) {
if (is.AFNI.1D.string(t)) return(eval.AFNI.1D.string(t))
else if (is.AFNI.R.string(t)) return(eval.AFNI.R.string(t))
else return(NULL)
}
date.stamp <- function (fancy=FALSE) {
if (fancy) {
return(gsub(' ','_',date()))
} else {
return(format(Sys.time(), "%d%m%y-%H%M%S"))
}
}
parse.AFNI.name.selectors <- function(filename,verb=0) {
n <- list()
n$brsel<- NULL
n$rosel<- NULL
n$rasel<- NULL
n$insel<- NULL
selecs <- strsplit(filename,"\\[|\\{|<|#")[[1]]
n$name <- selecs[1]
for (ss in selecs[2:length(selecs)]) {
if (length(grep("]",ss))) {
n$brsel <- strsplit(ss,"\\]")[[1]][1]
} else if (length(grep("}",ss))) {
n$rosel <- strsplit(ss,"\\}")[[1]][1]
} else if (length(grep(">",ss))) {
n$rasel <- strsplit(ss,">")[[1]][1]
}
}
selecs <- strsplit(filename,"#")[[1]]
if (length(selecs) > 1) {
n$insel <- selecs[2]
}
return(n)
}
parse.AFNI.name <- function(filename, verb = 0) {
if (filename == '-self_test') { #Secret testing flag
note.AFNI('Function running in test mode')
show.AFNI.name(parse.AFNI.name('DePath/hello.DePrefix', verb))
show.AFNI.name(parse.AFNI.name('DePath/DePrefix+acpc', verb))
show.AFNI.name(parse.AFNI.name('DePath/DePrefix+acpc.', verb))
show.AFNI.name(parse.AFNI.name('DePath/DePrefix+acpc.HEAD', verb))
show.AFNI.name(parse.AFNI.name('DePath/DePrefix+acpc.BRIK.gz', verb))
show.AFNI.name(parse.AFNI.name('DePath/DePrefix+acpc.HEAD[23]', verb))
show.AFNI.name(
parse.AFNI.name('DePath/DePrefix+acpc.HEAD[DeLabel]{DeRow}', verb))
show.AFNI.name(
parse.AFNI.name('DePath/DePrefix+acpc[DeLabel]{DeRow}', verb))
show.AFNI.name(
parse.AFNI.name('DePath/DePrefix+acpc.[DeLabel]{DeRow}', verb))
return(NULL)
}
an <- list()
an$view <- NULL
an$pprefix <- NULL
an$brsel <- NULL
an$rosel <- NULL
an$rasel <- NULL
an$insel <- NULL
an$type <- NULL
an$path <- NULL
an$orig_name <- filename
an$file <- NULL
if (verb) { cat ('Parsing >>',filename,'<<\n', sep=''); }
if (!is.character(filename)) {
warning(paste('filename >>',
filename, '<< not a character string\n', sep=''),
immediate. = TRUE)
traceback()
return(NULL)
}
#Deal with special names:
if (length(grep("^1D:.*$",filename))) {
an$type = '1Ds'
return(an)
} else if (length(grep("^R:.*$",filename))) {
an$type = 'Rs'
return(an)
}
#Deal with selectors
n <- parse.AFNI.name.selectors(filename, verb)
filename <- n$name
an$file <- n$name
an$brsel <- n$brsel
an$rosel <- n$rosel
an$rasel <- n$rasel
an$insel <- n$insel
#Remove last dot if there
filename <- sub('\\.$','',filename)
#NIFTI?
n <- strip.extension(filename, c('.nii', '.nii.gz'), verb)
if (n$ext != '') {
an$ext <- n$ext
an$type <- 'NIFTI'
an$pprefix <- n$name_noext
} else {
#remove other extensions
n <- strip.extension(filename, c('.HEAD','.BRIK','.BRIK.gz',
'.BRIK.bz2','.BRIK.Z',
'.1D', '.1D.dset',
'.niml.dset',
'.' ),
verb)
if (n$ext == '.1D' || n$ext == '.1D.dset') {
an$type <- '1D'
} else if (n$ext == '.niml.dset') {
an$type <- 'NIML'
} else {
an$type <- 'BRIK'
}
if (n$ext == '.') {
n$ext <- ''
}
an$ext <- n$ext
filename <- n$name_noext
n <- strip.extension(filename, c('+orig','+tlrc','+acpc'), verb)
if (n$ext != '') {
an$view <- n$ext
} else {
an$view <- NA
}
an$pprefix <- n$name_noext
}
#a prefix with no path
an$prefix <- basename(an$pprefix)
#and the path
an$path <- dirname(an$orig_name)
if (verb > 2) {
note.AFNI("Browser not active")
# browser()
}
if ( an$type != '1D' && (
!is.null(an$brsel) || !is.null(an$rosel) ||
!is.null(an$rasel) || !is.null(an$insel))) {
#Remove trailing quote if any
an$prefix <- gsub("'$", '', an$prefix)
an$prefix <- gsub('"$', '', an$prefix)
an$pprefix <- gsub("'$",'', an$pprefix)
an$pprefix <- gsub('"$','', an$pprefix)
}
if ( an$type != 'BRIK' ) {
#Put the extension back on
an$pprefix <- paste(an$pprefix,an$ext, sep='')
an$prefix <- paste(an$prefix,an$ext, sep='')
}
return(an)
}
exists.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an)
ans <- 0
if (file.exists(head.AFNI.name(an))) ans <- ans + 1;
if (file.exists(brik.AFNI.name(an)) ||
file.exists(paste(brik.AFNI.name(an),'.gz', sep='')) ||
file.exists(paste(brik.AFNI.name(an),'.Z', sep=''))) ans <- ans + 2
return(ans)
}
#Read parameters for -dataTable options (e.g. 3dLME)
#If opts is just one string starting with '@' such as
#@PPP, the parameters are read from text file PPP and
#returned as if they were found on the command line
dataTable.AFNI.parse <- function(opts) {
if (is.null(opts) || length(opts) == 0) {
return(NULL)
}
if (length(opts) == 1 && length(grep('^@.*$',opts))) {
ff <- strsplit(opts, '')[[1]]
ff <- paste(ff[2:length(ff)], sep='', collapse='')
if (!file.exists(ff)) {
err.AFNI(paste("table file", ff, "does not exist"));
return(NULL)
}
#Can't read as table because users might have EOL \ from shell command
#opts <- scan(ff, what='character')
opts <- scan(ff, what='character', na.strings="") # GC 08/15/2014: modified to allow NA in dataTable
opts <- opts[grep('[^\\\\]',opts)]
}
return(opts);
}
used.AFNI.prefix <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
ans <- 0
ans <- exists.AFNI.name(an)
if (ans > 0) return(1)
if (exists.AFNI.name(sprintf('%s+orig', an$pprefix))) return(1)
if (exists.AFNI.name(sprintf('%s+tlrc', an$pprefix))) return(1)
if (exists.AFNI.name(sprintf('%s+acpc', an$pprefix))) return(1)
return(0)
}
pprefix.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
return(an$pprefix);
}
prefix.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
return(an$prefix);
}
view.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
return(an$view);
}
pv.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
return(paste(an$pprefix,an$view,sep=''));
}
head.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
if (an$type == 'BRIK' && !is.na(an$view)) {
return(paste(an$pprefix,an$view,".HEAD",sep=''));
} else {
return((an$orig_name));
}
}
brik.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
if (an$type == 'BRIK' && !is.na(an$view)) {
return(paste(an$pprefix,an$view,".BRIK",sep=''));
} else {
return((an$orig_name));
}
}
compressed.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
if (length(grep('\\.gz$', an$ext))) {
return('gz')
} else if (length(grep('\\.bz2$', an$ext))) {
return('bz2')
} else if (length(grep('\\.Z$', an$ext))) {
return('Z')
} else {
return('')
}
}
modify.AFNI.name <- function (name, what="append", val="_new", cwd=NULL) {
if (!is.loaded('R_SUMA_ParseModifyName')) {
err.AFNI("Missing R_io.so");
return(NULL);
}
an <- .Call("R_SUMA_ParseModifyName",
name = name,
what = what,
val = val,
cwd = cwd)
return(an)
}
AFNI.command.history <- function(ExecName=NULL, args=NULL, ohist=NULL) {
if (!is.loaded('R_SUMA_HistString')) {
return(NULL);
}
an <- .Call("R_SUMA_HistString", ExecName, args, ohist)
return(an)
}
uncompress.AFNI <- function(an, verb = 1) {
if (is.character(an)) an <- parse.AFNI.name(an);
ans <- 0
zz <- paste(brik.AFNI.name(an),'.Z', sep='');
if (file.exists(zz)) {
if (verb) {
cat ('Uncompressing', zz, '\n');
}
system(paste('uncompress', zz));
}
zz <- paste(brik.AFNI.name(an),'.gz', sep='');
if (file.exists(zz)) {
if (verb) {
cat ('gzip -d-ing', zz, '\n');
}
system(paste('gzip -d', zz));
}
zz <- paste(brik.AFNI.name(an),'.bz2', sep='');
if (file.exists(zz)) {
if (verb) {
cat ('bzip2 -d-ing', zz, '\n');
}
system(paste('bzip2 -d', zz));
}
return(an);
}
copy.AFNI.dset <- function(an, ancp, overwrite=FALSE, verb = 1) {
if (is.character(an)) an <- parse.AFNI.name(an);
if (is.character(ancp)) ancp <- parse.AFNI.name(ancp);
ans <- 0
if (!exists.AFNI.name(an)) {
err.AFNI(sprintf('Dset %s not found', an$file));
return(ans)
}
if (exists.AFNI.name(ancp) && !overwrite) {
err.AFNI(sprintf('New dset %s already exists', ancp$file));
return(ans)
}
opt <- ''
if (overwrite) opt <- sprintf('%s %s', opt, -overwrite);
com = sprintf('3dcopy %s %s', an$file, ancp$file);
system(com, ignore.stderr = TRUE, intern=TRUE)
if (!exists.AFNI.name(ancp)) {
err.AFNI(sprintf('New dset %s not created', ancp$file));
return(ans)
}
return(1)
}
show.AFNI.name <- function(an) {
cat ('\n',
'uname=', an$orig_name, '\n',
'file=', an$file,'\n',
'type =', an$type,'\n',
'pref.=', an$pprefix, '\n',
'view =', an$view, '\n',
'head =', head.AFNI.name(an), '\n',
'brik =', brik.AFNI.name(an), '\n',
'brsel=', an$brsel, '\n',
'rosel=', an$rosel, '\n',
'rasel=', an$rasel, '\n',
'insel=', an$insel, '\n',
'compr=', compressed.AFNI.name(an), '\n',
'exist=', exists.AFNI.name(an), '\n');
}
#------------------------------------------------------------------
# Functions to parse command-line arguments
#------------------------------------------------------------------
value.AFNI.args <- function(name, ops) {
ifnd <- which(name == names(ops));
if (length(ifnd)) {
vv <- vector(typeof(ops[ifnd][[1]]));
for (i in 1:length(ifnd)) {
vv <- c (vv,ops[i][[1]])
}
return(vv);
}
return(NULL);
}
show.AFNI.args <- function (ops, verb=0, adieu=FALSE, hstr='') {
if (is.null(ops)) {
cat ('NULL options\n');
} else {
cat (hstr,'Allowed Options:\n');
if (length(ops[['allowed_options']])) {
for (i in 1:length(ops[['allowed_options']])) {
cat (' ', ops[['allowed_options']][i], '\n');
}
} else {
cat ('whatever grinds your beans');
}
if (verb) {
cat (hstr, 'User options:\n');
for (i in 1:length(ops)) {
if ((names(ops)[i] != 'allowed_options')) {
cat (' ', names(ops)[i], ': ', ops[[i]],'\n');
}
}
}
}
if (adieu) exit.AFNI(0)
}
check.AFNI.args <- function ( ops, params = NULL, verb=0) {
if (!is.null(params) && !is.null(ops)) {
for (i in 1:length(ops)) {
if (verb) {
str(names(ops)[i])
str(params[names(ops)[i]][[1]])
cat('\nChecking on ', paste(ops[[i]],collapse=','),'\n');
}
ipar <- which(names(ops)[i] == names(params));
if (length(ipar)) {
pp <- params[ipar[1]][[1]]['count'][[1]];
opsvec <- ops[[i]];
if (length(pp) == 1) { #exact number
if (length(opsvec) != pp) {
#browser()
msg <- paste( 'Expecting ',pp, ' parameters for option "',
names(ops)[i], '".\n Have ',
length(opsvec), ' parameter(s) in string "',
paste(opsvec, collapse = ' '),
'" instead.', sep = '')
if (length(opsvec) > pp &&
length(grep('^-[a-z,A-Z]', opsvec[1+pp]))) {
msg <- paste( msg, '\n NOTE that ', opsvec[1+pp],
' in bad option above is not a recognized option.',
collapse = '', sep = '' );
}
err.AFNI(msg);
return(0);
}
} else if (length(pp) == 2) { #range
if (length(opsvec) < pp[1] || length(opsvec) > pp[2]) {
if (pp[2] == Inf) {
msg <- paste( 'Expecting more than ',pp[1],
' parameters for option "',
names(ops)[i], '".\n Have ',
length(opsvec), ' parameter(s) in string "',
paste(opsvec, collapse = ' '),
'" instead.', sep = '')
if (length(opsvec) > 0 &&
length(grep('^-[a-z,A-Z]', opsvec[1]))) {
msg <- paste( msg, '\n NOTE that ', opsvec[1],
' in bad option above is not a recognized option.',
collapse = '', sep = '' );
}
err.AFNI(msg);
} else {
msg <- paste( 'Expecting ',pp[1], ' to ', pp[2],
' parameters for option "',
names(ops)[i], '".\n Have ',
length(opsvec), ' parameter(s) in string "',
paste(opsvec, collapse = ' '),
'" instead.', sep = '');
if (length(opsvec) > pp[2] &&
length(grep('^-[a-z,A-Z]', opsvec[1+pp[2]]))) {
msg <- paste( msg, '\n NOTE that ', opsvec[1+pp[2]],
' in bad option above is not a recognized option.',
collapse = '', sep = '' );
}
err.AFNI(msg);
}
return(0);
}
} else {
warning(paste( 'I do not know what to do here'),
immediate. = TRUE);
return(0);
}
} else {
#ok
}
}
}
return(1); #all ok
}
apl <- function ( n = 0, d=NA, h=NULL, dup=FALSE ) {
return(list('count'=n, 'default'=d, help=h, duplicate_ok=dup));
}
load.debug.AFNI.args <- function ( fnm = NULL) {
if (is.null(fnm)) {
fnm <- system('ls *.dbg.AFNI.args .*.dbg.AFNI.args', ignore.stderr = TRUE, intern=TRUE)
}
if (length(fnm) > 1) {
err.AFNI(paste("More than one argument file found:",
paste(fnm, collapse=' ', sep=' '), ".\n",
"Try one of:\n",
paste("load.debug.AFNI.args('",fnm,"')",
collapse='\n', sep=''),
sep = ''));
return(FALSE);
}else if (length(fnm) == 0) {
err.AFNI("No files to load");
return(FASLE);
}
load(fnm)
ok_auto <- FALSE
if (exists('rfile')) {
ss <- sprintf("\n Start interactive code with:\n\nsource('%s')\n",
rfile);
ok_auto <- TRUE
} else {
ss <-
sprintf("\n Start interactive code with:\n\nsource('yourRfile')\n");
}
if (exists('args')) {
note.AFNI(sprintf("Setting .DBG_args from args in %s%s", fnm, ss));
.DBG_args <<- args
} else {
err.AFNI(sprintf("Variable args not in %s", fnm));
return(FALSE)
}
if (ok_auto) {
if (prompt.AFNI("Execute source command?", c('y','n')) == 1) {
source(rfile)
}
}
return(TRUE)
}
parse.AFNI.args <- function ( args, params = NULL,
other_ok=TRUE,
verb = 0) {
#for (i in 1:length(args)) {
# cat (i, args[[i]],'\n');
#}
if (is.null(args)) return(NULL)
if (!is.null(params)) {
allowed_options <- sort(names(params));
duplicate_okvec <- vector('character');
for (i in 1:1:length(params)) {
pl <- params[i][[1]];
if (pl['duplicate_ok'][[1]]) {
duplicate_okvec <- c(duplicate_okvec, names(params)[i])
}
}
} else {
allowed_options <- vector('character');
}
#Add global args like -overwrite
allowed_options <- c(allowed_options, "-overwrite")
#find locations of -*
ii <- grep ('^-.*', args);
iflg <- vector('numeric')
if (length(ii) > 0) {
for (i in 1:length(ii)) {
if (!is.num.string(args[[ii[i]]])) {
if (!length(allowed_options)) {
iflg <- append(iflg, ii[i]);
} else { #Make sure it is an acceptable name
if (length(which(args[[ii[i]]] == allowed_options))) {
iflg <- append(iflg, ii[i]);
}
}
}
}
}
if (verb) note.AFNI(paste(args[iflg]))
ops = list()
used <- vector('logical', length(args));
if (length(iflg)) {
iflg <- append(iflg,length(args)+1)
#store results
nm <- vector('character');
for (i in 1:(length(iflg)-1)) {
if (0) { # Why remove the -?, makes things inconsistent elsewhere
#newnm <- strsplit(args[[iflg[i]]],'-')[[1]][2]
} else newnm <- args[[iflg[i]]]
if (length(nm) && length(which(newnm == nm)) &&
(newnm != '-gltLabel') && (newnm != '-gltCode') && # 10/18/2012 GC: added this line for 3dMVM
(newnm != '-glfLabel') && (newnm != '-glfCode') && # 12/22/2014 GC: added this line for 3dMVM
(!length(duplicate_okvec) ||
length(which(iflg[i] == duplicate_okvec))) ){
warning(paste('option ', newnm, 'already specified.\n'),
immediate. = TRUE);
show.AFNI.args(ops)
return(NULL);
}
#nm <- append(nm, newnm)
used[iflg[i]] <- TRUE;
istrt = iflg[i]+1;
pp <- vector('character');
if (istrt <= length(args) && istrt != iflg[i+1]) {
iend <- max(c(iflg[i+1]-1, istrt))
for (ii in istrt:iend) {
pp <- append(pp, args[[ii]]);
used[ii] <- TRUE;
}
}
#create a cleaned up string
pp <- paste(pp, collapse = ' ')
if ( length(grep('^".*"$',pp)) || #Quoted string, do not split
length(grep("^'.*'$",pp)) || #Quoted string, do not split
length(grep("^1D:.*$",pp)) || #1D: string, do not split
length(grep("^R:.*$",pp)) ) {
} else {
if (verb) {
note.AFNI(sprintf("Splitting >>>%s<<<", pp))
}
pp <- strsplit(clean.args.string(pp), ' ')
}
if((newnm != '-gltLabel') && (newnm != '-glfLabel') && (newnm != '-gltCode') && (newnm != '-glfCode')) { # 12/22/2014 GC: added this line for 3dMVM
ops <- c(ops, (pp))
names(ops)[length(ops)] <- newnm
} else if(length(which(newnm == nm)))
ops[[newnm]] <- c(ops[[newnm]], pp) else {
ops <- c(ops, list(pp))
names(ops)[length(ops)] <- newnm
}
nm <- append(nm, newnm)
}
}
#cleanup
if (length(ops)) {
for (i in 1:length(ops)) {
#ops[[i]] <- clean.args.string(ops[[i]])
}
}
#numeric changes
if (length(ops))
for (i in 1:length(ops))
if(!is.list(ops[[i]])) if (is.num.string(ops[[i]]))
ops[[i]] <- as.numeric(ops[[i]])
#defaults
pp <- c(args[used == FALSE])
ops <- c (ops, list("other"=pp));
#add allowed options
ops <- c (ops, list("allowed_options"=allowed_options));
if (!other_ok) {
if (length(ops[['other']])) {
err.AFNI(paste('Illegal parameters on command line:\n',
' ', ops['other'],
'\nTry -allowed_options, or -help for details\n',
'\n'));
exit.AFNI(1);
}
}
#check
if (!check.AFNI.args(ops, params)) {
return(NULL);
} else {
return(ops);
}
}
AFNI.new.options.list <- function(history = '', parsed_args = NULL) {
lop <- list (com_history = history);
#Look for defaults
lop$overwrite <- FALSE
for (i in 1:length(parsed_args)) {
opname <- strsplit(names(parsed_args)[i],'^-')[[1]];
opname <- opname[length(opname)];
switch(opname,
overwrite = lop$overwrite <- TRUE )
}
return(lop)
}
#------------------------------------------------------------------
# Some utilities
#------------------------------------------------------------------
#history.AFNI (is almost identical to history() but without
#the interactive mode and grep
history.AFNI <- function (max.show = 25, reverse = FALSE, pattern, ...)
{
file1 <- tempfile("Rrawhist")
savehistory(file1)
rawhist <- readLines(file1)
unlink(file1)
if (!missing(pattern))
rawhist <- unique(grep(pattern, rawhist, value = TRUE,
...))
nlines <- length(rawhist)
if (nlines) {
inds <- max(1, nlines - max.show):nlines
if (reverse)
inds <- rev(inds)
}
else inds <- integer(0L)
file2 <- tempfile("hist")
writeLines(rawhist[inds], file2)
#file.show(file2, title = "R History", delete.file = TRUE)
mm <- readLines(file2)
file.remove(file2)
cat(mm, sep='\n')
invisible(mm)
}
hgrep <- function (pattern=NULL){
if (is.null(pattern)) {
history.AFNI(max.show=200)
} else {
history.AFNI(max.show=Inf, pattern=pattern)
}
}
hsgrep <- function (pattern='source'){
hgrep(pattern)
}
#Report objects using the most memory
R.bit.version <- function () {
if (.Machine$sizeof.pointer == 4) return(32)
else if (.Machine$sizeof.pointer == 8) return(64)
else return(0)
}
memory.hogs <- function (n=10, top_frac=0.9, test=FALSE, msg=NULL) {
if (test) {
toy <- function() { g <- array(0,c(128,128,128,10));
memory.hogs(test=FALSE) }
moy <- function() {k <- array(0,c(128, 128, 128, 3)); toy() }
moy()
return(NULL)
}
us <- 1024*1024
nfr <- sys.nframe()
z<- NULL
pp <- 0
while (nfr > 0) {
nfr <- nfr -1
pp <- pp+1
envir <- parent.frame(pp)
envfunc <- (as.character(sys.call(-(pp))))
if (!length(envfunc)) envfunc <- 'GLOBAL'
zs <- sapply(ls(envir), function(x)
object.size(get(x, envir = envir)))
names(zs) <- paste(names(zs)," (in ",envfunc,")")
z <- c(z, zs)
}
#Stop at the top_frac
if (top_frac > 0.0) {
m <- rev(sort(z))[1:n]
mtot <- sum(m)
ms <- m[1]
nf <- 1
while (ms /mtot < top_frac && nf <= n) {
nf <- nf+1
ms <- ms+m[nf]
}
}
m <- as.matrix(rev(sort(z))[1:nf]/us)
colnames(m) <- c("Size(Mb)")
m <- rbind(m,' total' = sum(m));
m <- rbind(m,' Grand total' = sum(z)/us);
cat( sprintf("Memory hogs check from function: %s\nR-%d bit version\n",
as.character(paste(sys.call(-1), collapse=' ')),
R.bit.version()) );
if (!is.null(msg)) cat(msg)
print(m, digits=2)
invisible(m)
}
newid.AFNI <-function(len=26) {
return(paste("XYZ_",
paste(sample(c(rep(0:9,each=5),LETTERS, letters),
len-4, replace=TRUE), collapse='', sep=''), sep='',collapse=''))
}
sys.AFNI <- function(com=NULL, fout=NULL, ferr=NULL, echo=FALSE,
CoS=errex.AFNI ) {
ss <- list(stat=1, err='', out='')
if (is.null(com)) return(ss)
rmfout <- FALSE
if (is.null(fout)) {
fout <- sprintf('/tmp/fout.%s', newid.AFNI())
rmfout <- TRUE
}
rmferr <- FALSE
if (is.null(ferr)) {
ferr <- sprintf('/tmp/ferr.%s', newid.AFNI())
rmferr <- TRUE
}
if (echo) note.AFNI(paste("Executing shell command",com))
ss$stat<-try(system(paste(com,' 1>', fout,' 2>', ferr, collapse='')))
ss$out <- readLines(fout, n=-1, ok=TRUE, warn=FALSE, encoding='unknown')
if (rmfout) unlink(fout);
ss$err <- readLines(ferr, n=-1, ok=TRUE, warn=FALSE, encoding='unknown')
if (rmfout) unlink(ferr);
if (ss$stat && !is.null(CoS)) CoS(paste("Error status executing:\n",com))
invisible(ss)
}
who.called.me <- function (quiet_inquisitor=FALSE, trim = 0) {
mm <- (as.list(sys.calls()))
#str(mm)
N_mm <- length(mm)
callstr <- NULL
if (quiet_inquisitor) skp <- 2
else skp <- 1
callstr <- ''
for (i in (N_mm-skp):1 ) {
caller <- as.character(mm[i])
if (length(caller) == 0) {
caller <- 'R_prompt'
}
if (trim==-1) { #function only
caller <- strsplit(caller[1],'(', fixed=TRUE)[[1]][1]
} else if (trim>0) {
fun <- strsplit(caller[1],'(', fixed=TRUE)[[1]][1]
par <- strsplit(caller[1],'(', fixed=TRUE)[[1]][2]
par <- strsplit(par,'')[[1]]
n_char <- min(length(par),trim)
if (n_char < length(par)) ell <- '...'
else {
ell <- ''
#and remove last parentheses
par <- sub(')$','',par)
}
if (is.na(par[1])) {
caller <- paste(fun, collapse='', sep='')
} else {
caller <- paste(fun, '(', paste(par[1:n_char],collapse='', sep=''),
ell, ')', collapse='', sep='')
}
}
spc = ' '
if (i == N_mm-skp) {
callstr <- paste(callstr, caller[1])
} else {
spc <- paste(spc, ' ', sep='')
callstr <-paste(callstr, '\n',spc, '-->',
caller,
sep= '')
}
}
#if (BATCH_MODE) caller <- as.character(sys.call(-1))
#else caller <- as.character(sys.call(-2))
return(callstr)
}
#return 1 if all strings in vector ss can be changed to numbers
is.num.string <- function(ss) {
if (is.null(ss) || !length(ss) || ss == '' ||
is.null(tryCatch(as.numeric(ss),
warning=function(ex) {}))) {
return(0);
} else {
return(1);
}
}
clean.args.string <- function(ss) {
if (is.list(ss) || length(ss) > 1) {
warning(paste('Function only works on single strings',
str(ss),'\n', sep=''),
immediate.=TRUE);
return(NULL);
}
#remove trailing whites
ss <- sub('^[[:space:]]*','',ss);
ss <- sub('[[:space:]]*$','',ss);
#remove multiple whites
ss <- gsub('[[:space:]]+',' ',ss);
#treat = nicely
ss <- gsub('[[:space:]]*=[[:space:]]*','=',ss)
return(ss)
}
deblank.string <- function(s, start=TRUE, end=TRUE, middle=FALSE) {
if (end) {
s = sub('[[:space:]]+$','',s);
}
if (start) {
s = sub('^[[:space:]]+','',s);
}
if (middle) {
s = gsub('[[:space:]]+',' ',s);
}
return(s);
}
pad.string.lines <- function(s, pre=' ', post=NULL) {
if (!is.null(pre)) {
s = sub('^',pre,s);
s = gsub('\n',sprintf('\n%s',pre),s);
}
if (!is.null(post)) {
s = sub('$',post,s);
s = gsub('\n',sprintf('%s\n',post),s);
}
return(s);
}
trim.string <- function (s, nchar=32, left=TRUE, strim='...')
{
ss <- strsplit(s,'')[[1]]
if (length(ss)>nchar) {
#try deblanking
s <- deblank.string(s)
ss <- strsplit(s,'')[[1]]
nc <- length(ss)
if (nc>nchar) {
#browser()
nstrim = length(strsplit(strim,'')[[1]])
if (left) {
ns <- nc - nchar - nstrim
if (ns > nstrim) {
ss <- ss[ns:nc]
s<-paste(strim,paste(ss,collapse=''), sep='')
}
return(s)
}else {
ns <- nchar - nstrim
ss <- ss[1:ns]
s<-paste(paste(ss,collapse=''), strim, sep='')
}
} else return(s)
} else return(s)
}
as.num.vec <- function(ss, addcount=TRUE, sepstr='.', reset=FALSE) {
if (is.list(ss) || length(ss) > 1) {
warning(paste('Function only works on single strings',
str(ss),'\n', sep=''),
immediate.=TRUE);
return(NULL);
}
ss <- clean.args.string(ss)
dd <- strsplit(ss,' ')[[1]];
nn <- vector('numeric');
ww <- vector('character');
lastname <- '.v'
valnum <- 0
for (ii in 1:length(dd)) {
vv <- strsplit(dd[ii],'=')[[1]];
if (length(vv) > 1) {
valnum <- valnum+1
ll <- vv[1]
vv <- as.numeric(vv[length(vv)]);
if (is.na(vv)) { return(NULL); }
lastname <- ll
} else {
valnum <- valnum+1
wrn <- getOption('warn'); options(warn=-1);
vv <- as.numeric(vv[1]); options(warn=wrn);
if (is.na(vv)) { return(NULL); }
if (addcount) {
sfnd <- paste(lastname, sepstr,'[[:digit:]]*$', sep='', collapse='')
if (!reset) {
ifnd <- grep(sfnd,ww);
} else {
ifnd <- grep(sfnd,ww[length(ww)]);
if (length(ifnd)) {
ifnd <- length(ww);
} else {
valnum <- 1
}
}
if (length(ifnd)) {
lastval <- strsplit(ww[ifnd[length(ifnd)]],
paste(lastname, sepstr,sep=''))[[1]];
if (lastval[length(lastval)] == '') {
valnum <- 1
} else {
valnum <- as.numeric(lastval[length(lastval)]) + 1
}
}
ll <- paste(lastname,sepstr, as.numeric(valnum), sep='')
} else {
ll <- paste(lastname, sep='')
}
}
nn <- c(nn,vv)
ww <- c(ww,ll)
}
names(nn) <- ww
return(nn)
}
as.char.vec <- function(ss) {
if (is.list(ss) || length(ss) > 1) {
warning(paste('Function only works on single strings',
str(ss),'\n', sep=''),
immediate.=TRUE);
return(NULL);
}
ss <- clean.args.string(ss)
dd <- strsplit(ss,' ')[[1]];
nn <- vector('character');
ww <- vector('character');
for (ii in 1:length(dd)) {
vv <- strsplit(dd[ii],'=')[[1]];
if (length(vv) > 1) ll <- vv[1]
else ll <- paste('v',as.character(vv[1]), sep='')
vv <- as.character(vv[length(vv)]);
if (is.na(vv)) { return(NULL); }
nn <- c(nn,vv)
ww <- c(ww,ll)
}
names(nn) <- ww
return(nn)
}
#------------------------------------------------------------------
# Functions to read 1D and other tables
#------------------------------------------------------------------
expand_1D_string <- function (t) {
vvf = vector(length = 0, mode="numeric")
#replace .. with : and split at ,
s = strsplit(sub("..",":",t, fixed=TRUE), ",")[[1]]
#Now loop and form vector of components
for (si in s) {
#cat ("working ", si, "\n")
if (length(grep(":",si))) {
vv <- eval(parse(text=si))
} else if (length(grep("@",si))) {
ssi = as.numeric(strsplit(si,"@")[[1]])
#cat ("ssi = ", ssi, "\n")
vv <- rep(ssi[2], ssi[1])
} else {
vv <- as.numeric(si)
}
#cat(si," = ",vv, "\n")
vvf <- c(vvf, vv)
#cat("vvnow = ",vvf, "\n")
}
return(vvf)
}
is.wholenumber.AFNI <- function(x, tol = .Machine$double.eps^0.5, all=TRUE) {
if (is.null(x)) return(FALSE)
if (!all) {
return(abs(x - round(x)) < tol)
} else {
return(prod(abs(x - round(x)) < tol)==1)
}
return(FALSE)
}
r.NI_new_element <- function (name, dat=NULL, atlist=NULL, tp=NULL){
nel <- list (name=name, atlist=atlist, dat=NULL)
if (is.null(tp)) {
if (is.character(dat)) {
tp <- "String"
} else if (is.numeric(dat)) {
if (is.wholenumber.AFNI(dat)) tp<-"int"
else tp <- "float"
}
}
if (!is.null(dat)) {
if (is.vector(dat)) {
nel <- r.NI_set_attribute(nel, "ni_dimen", length(dat))
nel <- r.NI_set_attribute(nel, "ni_type", tp)
nel$dat <- matrix(as.character(dat), length(dat),1)
} else if (is.matrix(dat)){
nel <- r.NI_set_attribute(nel, "ni_dimen", paste(dim(dat),collapse=','))
nel <- r.NI_set_attribute(nel, "ni_type",
paste(rep(tp,dim(dat)[2]), collapse=','))
nel$dat <- matrix(as.character(dat), dim(dat)[1], dim(dat)[2])
} else {
err.AFNI("Bad dat");
}
if (tp == "String") {
nel$dat <- apply(nel$dat, c(1,2), paste)
}
}
return(nel)
}
r.NI_get_attribute <- function (nel,name, brsel=NULL,
colwise=FALSE, is1Dstr=FALSE,
sep=' ; ', num=FALSE) {
ffs <- NULL
for (i in 1:length(nel$atlist)) {
if (!is.na(nel$atlist[[i]]$lhs) &&
nel$atlist[[i]]$lhs == name) {
ffs <- nel$atlist[[i]]$rhs
break
}
}
if (is.null(ffs)) return(NULL)
#ffs <- gsub("[;\"]","", ffs)
if (colwise) { #a per column deal, process it
#remove bad chars and split into one for each column
if (is1Dstr) {
ffs = expand_1D_string (ffs)
} else {
ffsv <- strsplit(ffs,sep)[[1]]
#some components are blanks to be killed
ffs <- ffsv[which(nchar(ffsv)!=0)]
}
if (!is.null(brsel)) {
if (max(brsel+1) > length(ffs)) {
err.AFNI(paste("Have ", length(ffs),
"attribute elements in ", paste(ffs,collapse=' '),
".\nBrick selection calls for max. of ",
max(brsel+1)," columns\n"
));
}
ffs <- ffs[brsel+1]
if (num) return(as.numeric(ffs))
else return(ffs)
} else {
#No selection, return all
if (num) return(as.numeric(ffs))
else return(ffs)
}
} else {
if (num) return(as.numeric(ffs))
else return(ffs)
}
return(NULL)
}
r.NI_set_attribute <- function (nel,name, val) {
ig <- -1
if (length(nel$atlist)) {
for (i in 1:length(nel$atlist)) {
if (!is.na(nel$atlist[[i]]$lhs) &&
nel$atlist[[i]]$lhs == name) {
ig<-i
break
}
}
}
if (ig > 0) {
nel$atlist[[ig]] <- list(lhs=deblank.string(name),
rhs=r.NI_dequotestring(val, db=TRUE))
} else {
nel$atlist <- c (nel$atlist,
list(list(lhs=deblank.string(name),
rhs=r.NI_dequotestring(val, db=TRUE))))
}
return(nel)
}
r.NI_dequotestring <- function(val, db=FALSE) {
if (db) val <- deblank.string(val)
val <- sub('^\"','', val)
val <- sub('\"$','', val)
val <- sub("^\'",'', val)
val <- sub("\'$",'', val)
return((val))
}
# A very basic parser, works only on simple elements
#of ascii headersharp niml files . Needs lots of work!
r.NI_read_element <- function (fname, HeadOnly = TRUE) {
fnp <- parse.AFNI.name(fname)
if (!file.exists(fnp$file)) {
return(NULL)
}
ff <- scan(fnp$file, what = 'character', sep = '\n', quiet=TRUE)
return(r.NI_read_str_element(ff, HeadOnly))
}
r.NI_read_str_element <- function (ff, HeadOnly = TRUE) {
nel <- list(atlist=NULL, dat=NULL)
#Remove #
ff <- gsub('^[[:space:]]*#[[:space:]]*', '',ff)
if (!length(grep('^<', ff))) {
return(NULL)
}
#strt markers
strv <- grep ('<', ff)
#stp markers
stpv <- grep ('>', ff)
#check
if (length(strv) != length(stpv)) {
err.AFNI("Have extra brackets");
return(NULL)
}
if (length(strv) < 1) {
err.AFNI("Have nothing");
return(NULL)
}
#If you have a group, jump to next element.
shead <- ff[strv[1]:stpv[1]]
if (length(grep('ni_form[[:space:]]*=[[:space:]]*\"ni_group\"', shead))
&& stpv[1]+1<length(ff)) {
#try again
return(r.NI_read_str_element(ff[stpv[1]+1:length(ff)], HeadOnly))
}
#Get only the first element, for now
for (i in 1:1) {
shead <- ff[strv[i]:stpv[i]]
#get the name
nel$name <- strsplit(shead[1],'<')[[1]][2]
nel$atlist <- vector()
#remove last > from shead
shead[length(shead)] <- sub('>$','',shead[length(shead)])
for (j in 2:length(shead)) {
attr <- strsplit(shead[j],'=')[[1]]
if (length(attr) == 1) {
nel <- r.NI_set_attribute(nel, attr[1], "")
} else if (length(attr) == 2) {
nel <- r.NI_set_attribute(nel, attr[1], attr[2])
} else if (length(attr) > 2) {
err.AFNI(paste("Parse error for ", shead[j]));
}
}
}
if (HeadOnly) return(nel)
#Get the data part, the lazy way
#All the content here is text. Not sure what to do with
#this next. Wait till context arises
tp <- r.NI_get_attribute(nel,"ni_type");
if (tp == 'String') {
for (i in (stpv[1]+1):(strv[2]-1)) {
if (is.null(nel$dat))
nel$dat <- ff[i]
else nel$dat <- c(nel$dat, ff[i])
}
} else {
for (i in (stpv[1]+1):(strv[2]-1)) {
#browser()
if (is.null(nel$dat))
nel$dat <- rbind(r.NI_dequotestring(
strsplit(deblank.string(ff[i]),split=' ')[[1]]))
else nel$dat <- rbind(nel$dat,
r.NI_dequotestring(strsplit(
deblank.string(ff[i]),split=' ')[[1]]))
}
}
return(nel)
}
r.NI_write_str_element<- function(nel, HeadOnly = TRUE) {
ff <- paste("<",nel$name,'\n', sep='')
for (i in 1:length(nel$atlist)) {
ff <- paste(ff,
' ',nel$atlist[[i]]$lhs,'="', nel$atlist[[i]]$rhs,'"\n', sep='')
}
ff <- paste(ff,'>');
if (!HeadOnly) {
if (!is.null(nel$dat)) {
ff <- paste(ff,'\n')
if (is.matrix(nel$dat)) {
for (j in 1:dim(nel$dat)[1]) {
ff <- paste(ff,paste(nel$dat[j,], collapse= ' ', sep=''),
'\n', sep='')
}
} else {
ff <- paste(ff,paste(nel$dat, collapse= '', sep=''),'\n',sep='')
}
}
}
ff <- paste(ff,'</',nel$name,'>', sep='')
return(ff)
}
is.NI.file <- function (fname, asc=TRUE, hs=TRUE) {
fnp <- parse.AFNI.name(fname)
if (!file.exists(fnp$file)) {
return(FALSE)
}
ff <- scan(fnp$file, what = 'character', nmax = 2, sep = '\n', quiet=TRUE)
if (asc && hs) {
if (!length(grep('^[[:space:]]*#[[:space:]]*<', ff))) {
return(FALSE)
} else {
return(TRUE)
}
}
return(FALSE)
}
apply.AFNI.matrix.header <- function (fname, mat,
brsel=NULL, rosel=NULL,rasel=NULL,
nheadmax = 10, checkNI=TRUE) {
attr(mat,'name') <- 'noname'
attr(mat,'FileName') <- fname
fnp <- parse.AFNI.name(fname)
if (!file.exists(fnp$file)) {
return(mat)
}
#Does this look 1D nimly?
if (checkNI && !is.NI.file(fnp$file, asc=TRUE, hs=TRUE)) {
return(mat)
}
nel <- r.NI_read_element(fnp$file, HeadOnly = TRUE)
# Common attributes
attr(mat,'name') <- nel$name
if (!is.null(labels <-
r.NI_get_attribute(nel, 'ColumnLabels', brsel, colwise=TRUE))){
if (length(labels) == ncol(mat)) colnames(mat) <- labels
}
if (nel$name == 'matrix') {
if (!is.null(colg <-
r.NI_get_attribute(nel, 'ColumnGroups', brsel,
colwise=TRUE, is1Dstr=TRUE))){
if (length(colg) == ncol(mat)) attr(mat,'ColumnGroups') <- colg
if (!is.null(labels)) {
llv <- vector(length=0,mode="numeric")
tcolg <- unique(colg)
for (i in tcolg) {
if (i>0) {
ll <- which(colg == i)
llv <- c(llv,ll[1])
}
}
tnames <-paste(strsplit(labels[llv],"#0"))
attr(mat,'TaskNames') <- tnames
}
}
if (!is.null(TR <- r.NI_get_attribute(nel, 'RowTR'))){
attr(mat, 'TR') <- as.double(TR)
}
} else if (nel$name == 'DICE') {
} else if (nel$name == '3dhistog' || nel$name == 'seg_histogram') {
if (!is.null(bw <- r.NI_get_attribute(nel, 'BinWidth'))){
attr(mat, 'BinWidth') <- as.double(bw)
}
if (!is.null(bw <- r.NI_get_attribute(nel, 'window'))){
attr(mat, 'BinWidth') <- as.double(bw)
}
if (!is.null(bw <- r.NI_get_attribute(nel, 'xlabel'))){
attr(mat, 'xlabel') <- bw
}
} else {
if (0) { #No need to whine
warn.AFNI(paste(
"Don't know what to do with attribute of element ", nel$name));
}
}
return(mat)
}
read.AFNI.xmat <- function (xmatfile, nheadmax=10) {
if (!file.exists(xmatfile)) {
err.AFNI(paste("xmatfile ",xmatfile,"not found"));
return(NULL);
}
#Now load the whole deal, comments are skipped by defaults
ff <- as.matrix(read.table(xmatfile))
#Check for attributes
ff <- apply.AFNI.matrix.header(xmatfile, ff)
#Make matrix be ts
ff <- ts(ff, names=colnames(ff), deltat=attr(ff,"TR"))
return(ff)
}
xmat.base.index <- function (xmat, AFNIindex = FALSE) {
cg <- attr(xmat,'ColumnGroups')
an <- which(cg == -1)
if (AFNIindex) an <- an -1
return(an)
}
xmat.motion.index <- function (xmat, AFNIindex = FALSE) {
cg <- attr(xmat,'ColumnGroups')
an <- which(cg == 0)
if (AFNIindex) an <- an -1
return(an)
}
xmat.roni.index <- function (xmat, AFNIindex = FALSE) {
cg <- attr(xmat,'ColumnGroups')
an <- which(cg <= 0)
if (AFNIindex) an <- an -1
return(an)
}
xmat.alltasks.index <- function (xmat, AFNIindex = FALSE) {
cg <- attr(xmat,'ColumnGroups')
an <- which(cg > 0)
if (AFNIindex) an <- an -1
return(an)
}
xmat.select.indices <- function(selstr, xmat, AFNIindex = FALSE) {
#paste(selstr, collapse=" ") is used because at times selstr is a vector
sel <- strsplit(paste(selstr, collapse=" "),
" ")[[1]]
#str(sel)
#cat(sel, 'length:' , length(sel), '\n')
ilst = vector('integer');
if (length(sel)) {
for (isel in sel) {
#cat("processing: ",isel,'\n')
if (!inherits( e <- try(as.integer(eval(parse(text=isel))),
silent=TRUE),
"try-error")) {
ilst <- append(ilst,
e+1,
after=length(ilst))
} else {
#cat("processing as string: ",isel,'\n')
if (isel == 'ALL') {
ilst <- append(ilst, seq(1,ncol(xmat)))
} else if (isel == 'ALL_TASKS') {
ilst <- append(ilst,
xmat.alltasks.index(xmat, AFNIindex=FALSE))
} else if (isel == 'RONI') {
ilst <- append(ilst, xmat.roni.index(xmat, AFNIindex=FALSE))
} else if (isel == 'MOTION') {
ilst <- append(ilst, xmat.motion.index(xmat, AFNIindex=FALSE))
} else if (isel == 'BASE') {
ilst <- append(ilst, xmat.base.index(xmat, AFNIindex=FALSE))
} else {
if (0) { #Too lose
ilst <- append(ilst, grep(isel, colnames(xmat)),
after=length(ilst))
} else {
ilst <- append(ilst,
grep(sprintf('^%s',isel), colnames(xmat)),
after=length(ilst))
}
}
}
}
} else {
ilst = 0:(ncol(xmat)-1)
}
if (length(ilst) && AFNIindex) ilst <- ilst -1
return(ilst)
}
read.AFNI.matrix.test <- function(verb=1) {
cat ( 'Running read.AFNI.matrix in test mode\n',
'See read.AFNI.matrix.test for details' )
i <- 0
mm <- paste ( 'subj age weight height\n',
'joe 13 299 123 \n',
'jane 22 600 234 \n',
'jim 2 188 23\n',
sep='', collapse='');
fouts <- sprintf('___fout%02d.1D', i);
cat (mm,file=fouts);
comm <- 'read.AFNI.matrix(fouts, verb=verb)'
note.AFNI(sprintf("Working with fouts=%s:\n %s", fouts, comm))
print(eval(parse(text=comm)))
comm <- paste("read.AFNI.matrix(fouts, verb=verb,",
"userrownames=c('jim','jane'))",
sep='', collapse='')
note.AFNI(sprintf("Working with fouts=%s:\n %s", fouts, comm))
print(eval(parse(text=comm)))
comm <- paste("read.AFNI.matrix(fouts, verb=verb, ",
"userrownames=c('jim','jane'),",
"usercolnames=c('weight','age'))",
sep='', collapse='')
print(eval(parse(text=comm)))
i<- i+1
mm <- paste ( ' age weight height\n',
' 13 299 123 \n',
' 22 600 234 \n',
' 2 188 23\n',
sep='', collapse='');
fouts <- sprintf('___fout%02d.1D', i);
cat (mm,file=fouts);
note.AFNI(sprintf("Working with %s", fouts))
comm <- 'read.AFNI.matrix(fouts, verb=verb)'
note.AFNI(sprintf("Working with fouts=%s:\n %s", fouts, comm))
print(eval(parse(text=comm)))
comm <- paste("read.AFNI.matrix(fouts, verb=verb,",
"usercolnames=c('height','weight,','age'))",
sep='', collapse='')
note.AFNI(sprintf("Working with fouts=%s:\n %s", fouts, comm))
print(eval(parse(text=comm)))
i<- i+1
mm <- paste ( ' 13 299 123 \n',
' 22 600 234 \n',
' 2 188 23\n',
sep='', collapse='');
fouts <- sprintf('___fout%02d.1D', i);
cat (mm,file=fouts);
note.AFNI(sprintf("Working with %s", fouts))
print(read.AFNI.matrix(fouts, verb=verb))
print(read.AFNI.matrix(fouts, verb=verb, usercolnames=c('col02'),
userrownames=c('row03','row01','row01')))
i<- i+1
mm <- paste ( ' heft \n',
' 299 \n',
' 600 \n',
' 188 \n',
sep='', collapse='');
fouts <- sprintf('___fout%02d.1D', i);
cat (mm,file=fouts);
note.AFNI(sprintf("Working with %s", fouts))
print(read.AFNI.matrix(fouts, verb=verb))
print(read.AFNI.matrix(fouts, verb=verb, usercolnames=c('heft')))
i<- i+1
mm <- paste ( ' 299 \n',
' 600 \n',
' 188 \n',
sep='', collapse='');
fouts <- sprintf('___fout%02d.1D', i);
note.AFNI(sprintf("Working with %s", fouts))
cat (mm,file=fouts);
print(read.AFNI.matrix(fouts, verb=verb))
i <- i+1
mm <- paste ( 'subj age \n',
'joe 13 \n',
'jane 22 \n',
'jim 2 \n',
sep='', collapse='');
fouts <- sprintf('___fout%02d.1D', i);
cat (mm,file=fouts);
comm <- 'read.AFNI.matrix(fouts, verb=verb)'
note.AFNI(sprintf("Working with fouts=%s:\n %s", fouts, comm))
print(eval(parse(text=comm)))
comm <- paste("read.AFNI.matrix(fouts, verb=verb,",
"userrownames=c('jim','jane'))",
sep='', collapse='')
note.AFNI(sprintf("Working with fouts=%s:\n %s", fouts, comm))
print(eval(parse(text=comm)))
for (j in 0:i) {
system(sprintf('\\rm -f ___fout%02d.1D', i));
}
}
AFNI.matrix.input.help.string <- function() {
s<-
'You specify 1D input in a variety of ways:
FILE.1D : Name of 1D file
1D_EXPR: A 1D expression a la AFNI "1D: 1,2,5"
R_EXPR: An R expression, "R: rep(seq(1,3),2)"
'
return(s)
}
read.AFNI.matrix <- function (fname,
usercolnames=NULL,
userrownames=NULL,
checkNI = TRUE, verb = 0) {
if (length(fname)>1) {
err.AFNI(paste("Cannot handle more than one name.\nHave ",
paste(fname,collapse=' '), sep=''))
return(NULL);
}
if (fname == '-self_test') {
read.AFNI.matrix.test()
return(NULL);
}
if (!is.null(mm <- eval.AFNI.string(fname))) {
#Might need to add some names someday...
attr(mm,'FileName') <- paste('STR<',fname,'>',sep='')
return(as.matrix(mm))
}
if (verb) cat(who.called.me())
if (is.character(fname)) {
fname <- parse.AFNI.name(fname)
fnameattr <- fname$orig_name
}else{
fnameattr <- 'NONAME'
}
if (fname$type == 'NIML') {
nmout <- sprintf('/tmp/fout.%s.1D.dset', newid.AFNI())
if (verb) warn.AFNI("Clumsy handling of NIML files via ConvertDset ...");
com <- paste(
'ConvertDset -overwrite -o_1D -prefix ', nmout, ' -i',
fname$orig_name);
ss <- sys.AFNI(com, echo = FALSE);
if (ss$stat == 1) {
err.AFNI(paste("Failed to get keys from labeltable",ltfile, ".\n",
"Command ",com,"Failed"));
return(NULL);
}
mm<-read.AFNI.matrix(nmout);
sys.AFNI(sprintf('\\rm -f %s',nmout));
return(mm);
} else {
#str(fname)
#fname$file
brk <- NULL
brk <- tryCatch({read.table(fname$file, colClasses='character')},
error=function(a){})
if (is.null(brk)) { #try as niml, just in case
if (verb) note.AFNI(paste("Attempting read as NIML."), callstr='');
fff <- r.NI_read_element(fname$file, FALSE)
if (!is.null(fff$dat)) {
brk <- as.data.frame(fff$dat, stringsAsFactors=FALSE)
checkNI = FALSE;
} else {
err.AFNI(paste("Failed to read matrix from ", fname$file,".\n"));
return(NULL);
}
}
}
if ( tolower(brk$V1[1]) == 'name' ||
tolower(brk$V1[1]) == 'subj' ||
tolower(brk$V1[1]) == '#file') {
subjCol <- brk$V1[2:dim(brk)[1]];
covNames <- paste(brk[1,2:dim(brk)[2]]);
if (dim(brk)[2]-1 == 1) {
covMatrix <- as.matrix(as.numeric(brk[2:dim(brk)[1],2]))
} else {
for (ii in 1:(dim(brk)[2]-1)) { #Add one column at a time
if (ii==1) {
covMatrix <- cbind(
as.numeric(brk[2:dim(brk)[1],2:dim(brk)[2]][[ii]]));
} else {
covMatrix <- cbind(covMatrix,
as.numeric(brk[2:dim(brk)[1],2:dim(brk)[2]][[ii]]));
}
}
}
} else {
flg <- tryCatch({as.numeric(brk$V1[1])}, warning=function(aa) {});
if (is.null(flg)) { #Just labels
covNames <- paste(brk[1,1:dim(brk)[2]]);
subjCol <- paste('row',sprintf('%02d',c(1:(dim(brk)[1]-1))), sep='')
istrt<- 2
}else { #completely naked
ccc <- trim.string(fname$prefix, nchar=10, left=FALSE, strim='..');
covNames <- paste(ccc,sprintf(' c%02d',c(1:dim(brk)[2])),sep='');
subjCol <- paste(ccc,sprintf(' r%02d',c(1:dim(brk)[1])), sep='')
istrt<- 1
}
if (dim(brk)[2] == 1) {
covMatrix <- as.matrix(as.numeric(brk[istrt:dim(brk)[1],1]))
} else {
for (ii in 1:(dim(brk)[2])) { #Add one column at a time
ccc <- tryCatch(
{as.numeric(brk[istrt:dim(brk)[1],1:dim(brk)[2]][[ii]])},
warning=function(aa) {} )
if (is.null(ccc)) {
warn.AFNI(paste("Failed to process column ",
ii-1, " in ", fnameattr,
". Using NA instead"))
ccc <- NA*vector(length=dim(brk)[1]-istrt+1)
}
if (ii==1) {
covMatrix <- cbind(ccc);
} else {
covMatrix <- cbind(covMatrix,ccc);
}
}
}
}
if (verb>2) {
note.AFNI("Browser here, not active");
#browser()
}
rownames(covMatrix) <- subjCol;
colnames(covMatrix) <- covNames;
#And now apply selectors
if (!is.null(fname$brsel)) {
sbsel <- eval.AFNI.1D.string(fname$brsel, nmax=dim(covMatrix)[2]-1)
if (min(sbsel) < 0 || max(sbsel)>=dim(covMatrix)[2]) {
err.AFNI(
sprintf('column selection outside possible range of <0,%d> in %s',
dim(covMatrix)[2]-1, fname$file));
return(NULL);
}
covMatrix <- covMatrix[,sbsel+1, drop=FALSE]
} else sbsel <- NULL
if (!is.null(fname$rosel)) {
rosel <- eval.AFNI.1D.string(fname$rosel, nmax=dim(covMatrix)[1]-1)
if (min(rosel) < 0 || max(rosel)>=dim(covMatrix)[1]) {
err.AFNI(
sprintf('row selection outside possible range of <0,%d> in %s',
dim(covMatrix)[1]-1, fname$file));
return(NULL);
}
covMatrix <- covMatrix[rosel+1,, drop=FALSE]
} else rosel <- NULL
if (!is.null(fname$rasel)) {
err.AFNI('Not ready to deal with range selection');
return(NULL);
} else rasel <- NULL
#Now, reorder per user*names
if (!is.null(userrownames)) {
dd <- userrownames[!(userrownames %in% subjCol)]
if (length(dd)) {
warning (paste('Row(s) "', paste(dd,collapse=' '),
'" do(es) not have an entry.\n'),
immediate.=TRUE);
return(NULL);
}
for (ii in 1:1:length(userrownames)) {
if (ii==1) {
mm <- rbind(covMatrix[userrownames[ii],]);
} else {
mm <- rbind(mm, covMatrix[userrownames[ii],]);
}
}
rownames(mm) <- userrownames
if(length(covNames)) colnames(mm) <- covNames
covMatrix <- mm
}
if (!is.null(usercolnames)) {
dd <- usercolnames[!(usercolnames %in% covNames)]
if (length(dd)) {
warning (paste('Column(s) "', paste(dd,collapse=' '),
'" do(es) not have an entry.\n'),
immediate.=TRUE);
return(NULL);
}
for (ii in 1:1:length(usercolnames)) {
if (ii==1) {
mm <- cbind(covMatrix[,usercolnames[ii]]);
} else {
mm <- cbind(mm, covMatrix[,usercolnames[ii]]);
}
}
colnames(mm) <- usercolnames
covMatrix <- mm
}
covMatrix <- apply.AFNI.matrix.header(fnameattr, covMatrix,
brsel=sbsel, rosel=rosel,
rasel=rasel,checkNI = checkNI)
return(covMatrix)
}
read.AFNI.labeltable <- function (ltfile=NULL,
def.grp.label = c('CSF','GM','WM','InSk','OutSk','Out')) {
if (!file.exists(ltfile)) {
err.AFNI(paste(ltfile, "does not exist"));
return(NULL)
}
if (is.null(lt <- r.NI_read_element(ltfile, HeadOnly=FALSE))) {
err.AFNI("Failed to read label table");
return(NULL)
}
if (lt$name != "VALUE_LABEL_DTABLE") {
err.AFNI("lt file not VALUE_LABEL_DTABLE")
return(NULL)
}
lt$labels <- lt$dat[,2]
lt$keys <- as.numeric(lt$dat[,1])
if (is.null(gl <- r.NI_get_attribute(lt, "grp.label", colwise=TRUE))) {
gl <- def.grp.label
}
if (!is.null(gl)) {
lt$grp.label <- gl
}
grps <- r.NI_get_attribute(lt, "grp.label", colwise=TRUE, num=TRUE)
if (is.null(grps)) {
for (lab in lt$labels) {
nn <- 0
ig <- 0
for (i in 1:length(gl)) {
if (length(grep(gl[i],lab))) {
if (ig == 0) {
nn <- 1
ig <- i
} else {
if ( length(gl[ig]) != length(lab) ) {
if (length(gl[i]) == length(lab) ) {
# A better match
nn <- 1
ig <- i
} else {
#Another partial
nn <- nn + 1
ig <- i
}
}
}
}
}
if (nn != 1) {
warn.AFNI(paste(
"Failed to find 1 group for ",lab, ". Found ",
nn ,"candidates. Setting grp to",
length(gl)+1))
ig <- length(gl)+1
}
grps <- c(grps, ig)
}
}
if (!is.null(grps)) {
lt$groups <- grps
}
return(lt)
}
write.AFNI.matrix <- function (m, fname='test.1D') {
if (class(m) == "AFNI_c_dataset") {
m <- m$brk
}
if (is.vector(m)) {
write(m, fname, ncolumns=1)
} else if (is.matrix(m)) { #Need to transpose matrix
write(t(m), fname, ncolumns=dim(m)[2])
} else {
err.AFNI("Don't know how to write this thing");
}
return(fname)
}
dset.nonzero.indices <- function (dset) {
nn <- dset.ndim(dset)
if (nn == 3 || nn==1) {
return(which(rowSums(abs(dset$brk), dims=nn)!=0,arr.ind=TRUE))
} else {
err.AFNI("Bad dim");
return(NULL);
}
}
dset.zero.indices <- function (dset) {
nn <- dset.ndim(dset)
if (nn == 3 || nn==1) {
return(which(rowSums(abs(dset$brk), dims=nn)==0,arr.ind=TRUE))
} else {
err.AFNI("Bad dim");
return(NULL);
}
}
dset.nvox <- function (dset) {
nvox = -1
if (length(dd <- dim(dset$brk)) <= 2) {
nvox = dd[1]
} else if (length(dd) > 2) {
nvox = prod(dd[1:3])
}
return(nvox)
}
dset.nval <- function (dset) {
nval = -1
if (length(dd <- dim(dset$brk)) == 1) {
nval = 1
} else if (length(dd) == 2) {
nval = dd[2]
} else if (length(dd) == 3) {
nval = 1
} else if (length(dd) == 4) {
nval = dd[4]
}
return(nval)
}
dset.nx <- function(dset) {
return(dset$dim[1])
}
dset.ny <- function(dset) {
return(dset$dim[2])
}
dset.nz <- function(dset) {
return(dset$dim[3])
}
dset.nslice <- function(dset) {
return(dset$dim[3])
}
dset.ndim <- function (dset) {
ndims <- -1
if (length(dd <- dim(dset$brk)) == 1 ||
length(dd) == 2) {
ndims <- 1
} else if (length(dd) == 3 || length(dd) == 4) {
ndims <- 3
}
return(ndims)
}
get.c.AFNI.attribute<- function (atlist, attribute, asis = FALSE, verb=0) {
if (!is.null(atlist$NI_head)) atlist <- atlist$NI_head
nel <- atlist[[attribute]]
if (!is.null(nel)) {
if (verb > 1) note.AFNI(c("Found", attribute, '\n'))
if (asis) { return(nel$dat) }
else {
tp <- r.NI_get_attribute(nel,"ni_type");
if (tp == 'String') {
dd <- sub("^ *\"",'', nel$dat)
dd <- sub(" *\"$",'', dd)
return(strsplit(dd,'~')[[1]])
} else if (tp == 'float') {
return(as.numeric(nel$dat));
} else if (tp == 'int') {
return(as.integer(nel$dat));
} else {
cat ("What of ", tp)
#browser()
nel
return(NULL);
}
}
}
if (verb) {
err.AFNI(paste("Attribute:",attribute,"Not found"))
show.dset.attr(atlist)
}
return(NULL);
}
set.c.AFNI.attribute<- function (atlist, attribute, val=NULL, verb=0, tp=NULL,
strsep="~") {
if (!is.null(atlist$NI_head)) atlist <- atlist$NI_head
nel <- atlist[[attribute]]
if (!is.null(nel)) {
if (verb) cat ("Found", attribute, '\n')
tp <- r.NI_get_attribute(nel,"ni_type");
if (tp == 'String') {
nel$dat <- paste("\"",paste(val, collapse=strsep), "\"", sep='')
nel <- r.NI_set_attribute(nel, "ni_dimen", 1);
} else if (tp == 'float' || tp == 'int') {
ff <- NULL
if (is.matrix(val)) {
nel <- r.NI_set_attribute(nel,
"ni_dimen",paste(dim(val),collapse=','))
nel$dat <- matrix(as.character(val), dim(val)[1], dim(val)[2])
} else if (is.vector(val)) {
nel <- r.NI_set_attribute(nel,"ni_dimen",paste(length(val)))
nel$dat <- matrix(as.character(val), length(val),1)
} else {
nel$dat <- val
}
} else {
cat ("What of ", tp)
#browser()
nel
}
} else {
nel <- r.NI_new_element("AFNI_atr", dat=val, atlist= NULL, tp=tp)
nel <- r.NI_set_attribute(nel, "atr_name", attribute)
atlist[[attribute]] <- nel
#Now that name and type are set, call again because
#you want the multi-col strings to be turned into one
#part with the '~' delimiters....
tp <- r.NI_get_attribute(nel,"ni_type");
if (tp == 'String') {
atlist <- set.c.AFNI.attribute(atlist, attribute, val=val,
verb=verb,strsep=strsep)
return(atlist);
}
}
atlist[[attribute]] <- nel
if (verb && is.null(nel)) {
err.AFNI(paste("Attribute:",attribute,"Not found."))
show.dset.attr(atlist)
}
return(atlist);
}
#See README.attributes's SCENE_DATA and niml_stat.c's distname[]
statsym.distold2new <- function(dist) {
if (dist == 'fico') return("Correl")
if (dist == 'fitt') return("Ttest")
if (dist == 'fift') return("Ftest")
if (dist == 'fizt') return("Zscore")
if (dist == 'fict') return("Chisq")
if (dist == 'fibt') return("Beta")
if (dist == 'fibn') return("Binom")
if (dist == 'figt') return("Gamma")
if (dist == 'fipt') return("Poisson")
return(dist)
}
statsym.list2code <- function(statsym) {
if (is.null(statsym) || length(statsym)==0) return("")
imx <- statsym[[1]]$sb
if(length(statsym)>1) for (i in 2:length(statsym)) {
if (imx < statsym[[i]]$sb) imx <- statsym[[i]]$sb
}
code <- NULL
for (i in 1:imx+1) {
code <- c(code,"none");
}
for (i in 1:length(statsym)) {
code[statsym[[i]]$sb+1] <-
paste(statsym.distold2new(statsym[[i]]$typ),"(",
paste(statsym[[i]]$par,collapse=','), ")",
sep = '');
}
return(code)
}
header.version <- function(hh, defmeth='UNKNOWN') {
if (defmeth == 'AUTO') {
if (have_R_io()) defmeth <- 'clib' else defmeth <- 'Rlib'
}
if (!is.null(hh) && length(hh)>0) {
if (is.null(names(hh[[1]]))) return('Rlib');
if (length(which("atlist" == names(hh[[1]]))) > 0) return('clib');
#Cannot tell, return default
return(defmeth)
}
return(defmeth)
}
#This function is to be used in different ways depending on whether you
#are setting, or retrieving attributes
#Example:
#rs <- read.c.AFNI('littleone.niml.dset')
#to retrieve attributes:
# dset.attr(rs,"IJK_TO_DICOM")
# or
# dset.attr(rs,"BRICK_LABS")
#To set attributes:
# rs <- dset.attr(rs,"BRICK_LABS",val=c('one', 'two', 'three'))
#To set an attribute to a default if none exists
# rs <- dset.attr(rs,"BRICK_LABS",default=c('set if none existed'))
#In query mode only, one can also pass the old '$header' list
#generated by the old read.AFNI
dset.attr <- function (dset, name=NULL, colwise=FALSE, num=FALSE,
val=NULL, default=NULL, tp=NULL) {
attr <- NULL
if (class(dset) != 'AFNI_R_dataset' && class(dset) != 'AFNI_c_dataset') {
if (header.version(dset) == 'Rlib') {
#OK to change nature of dset, only queries allowed for old style
dset <- list(header=dset);
}
}
if (!is.null(dset$header)) { #olde way
if (!is.null(val)) {
err.AFNI("Cannot set values with old header format");
return(attr)
}
iattr <- which(name == names(dset$header))
if (!length(iattr)) {
return(attr)
} else {
attr <- dset$header[[iattr]]
}
if (colwise) {
attr <- sub("^'",'', attr)
attr <- strsplit(attr,'~')[[1]]
}
if (num) attr <- as.numeric(attr)
} else {
if (!is.null(dset$NI_head)) hatr <- dset$NI_head
else hatr <- dset
if (is.null(name)) {
return(names(hatr))
}
if (is.null(val) && is.null(default)) { #Retrieval
if (name == 'dim') {
dd <- get.c.AFNI.attribute(hatr, "DATASET_DIMENSIONS")
ee <- get.c.AFNI.attribute(hatr, "DATASET_RANK")
return(c(dd[1], dd[2], dd[3], ee[2]))
} else if (name == 'orient') {
return(orcode.AFNI(get.c.AFNI.attribute(hatr, "ORIENT_SPECIFIC")))
} else if (name == 'origin') {
return(get.c.AFNI.attribute(hatr, "ORIGIN"))
} else if (name == 'delta') {
return(get.c.AFNI.attribute(hatr, "DELTA"))
} else if (name == 'hist') {
return(get.c.AFNI.attribute(hatr, "HISTORY_NOTE"))
} else if (name == "TR") {
dd <- get.c.AFNI.attribute(hatr, "TAXIS_NUMS")
ee <- get.c.AFNI.attribute(hatr, "TAXIS_FLOATS")
if (dd[3] == 77001) ee[2] = ee[2]/1000; #From msec to sec
return(ee[2])
} else {
return(get.c.AFNI.attribute(hatr, name))
}
} else if (!is.null(val)) { #Set the attribute
if (name == 'note') {
i <- get.c.AFNI.attribute(hatr, "NOTES_COUNT")
if (is.null(i)) i<-0
i<-i+1
hatr <- set.c.AFNI.attribute(hatr, "NOTES_COUNT", val = i)
name <- sprintf('NOTE_NUMBER_%03d', i)
hatr <- set.c.AFNI.attribute(hatr, name, val)
} else if (name == 'dim') {
if (is.null(val)) {
if (is.null(dset$brk)) {
err.AFNI("Cannot set dims from nothing");
return(NULL);
}
val <- dset.dimBRKarray(dset)
}
if (length(val) != 4) {
err.AFNI("Bad val for dim");
return(NULL);
}
hatr <- set.c.AFNI.attribute(hatr, "DATASET_DIMENSIONS",
val=c(val[1:3], 0, 0))
hatr <- set.c.AFNI.attribute(hatr, "DATASET_RANK",
val=c(3, val[4], rep(0,6)))
} else if (name == "statsym") {
hatr <- set.c.AFNI.attribute(hatr, "BRICK_STATSYM",
statsym.list2code(statsym=val), strsep=';')
} else if (name == "TR") {
dd <- get.c.AFNI.attribute(hatr, "TAXIS_NUMS");
if (is.null(dd)) {
rr <- get.c.AFNI.attribute(hatr, "DATASET_RANK")
if (is.null(rr)) {
err.AFNI("Cannot set TAXIS_NUMS with DATASET_RANK missing");
return(NULL);
}
dd <- c(rr[2], 0, 77002, -999, -999, -999, -999, -999);
}
ee <- get.c.AFNI.attribute(hatr, "TAXIS_FLOATS");
if (is.null(ee)) {
ee <- c(0, 0, 0, 0, 0, -999999.000, -999999.000, -999999.0000);
}
dd[2] <- 77002; #seconds
ee[2] <- val;
hatr <- set.c.AFNI.attribute(hatr, "TAXIS_NUMS", val=dd);
hatr <- set.c.AFNI.attribute(hatr, "TAXIS_FLOATS", val = ee);
} else{
hatr <- set.c.AFNI.attribute(hatr, name, val, tp=tp)
}
if (!is.null(dset$NI_head)) {
dset$NI_head <- hatr;
return(dset);
} else return(hatr);
} else if (!is.null(default)) {
if (is.null(get.c.AFNI.attribute(hatr,name))) {
hatr <- set.c.AFNI.attribute(hatr, name, val=default, tp=tp)
}
if (!is.null(dset$NI_head)) {
dset$NI_head <- hatr;
return(dset);
} else return(hatr);
}
}
return(attr)
}
show.dset.attr <- function (atlist, val=FALSE, sel=NULL ) {
if (!is.null(atlist$NI_head)) atlist <- atlist$NI_head
j = 1;
while (j<=length(atlist)) {
nm <- r.NI_get_attribute(atlist[[j]],"atr_name")
if (is.null(sel) || length(grep(sel,nm, ignore.case=TRUE))) {
cat(nm)
if (val) {
cat('=', atlist[[j]]$dat,'\n')
} else {
cat('\n');
}
}
j <- j + 1;
}
invisible(NULL);
}
dset.labels <- function (dset) {
return(dset.attr(dset, 'BRICK_LABS', colwise=TRUE))
}
dimBRKarray <- function(brk=NULL) {
d <- c(1,1,1,1)
if (is.array(brk)) dd <- dim(brk)
else if (is.vector(brk)) dd <- length(brk)
else {
note.AFNI("Don't know what to do with this brk")
return(NULL)
}
d[1:length(dd)]<-dd
return(d)
}
subBRKarray <- function(brk=NULL, sel=NULL) {
d <- dimBRKarray(brk)
if (is.null(d)) {
note.AFNI("Failed to get dims")
return(NULL)
}
if (is.null(sel) || d[4] == 1) {
v <- brk
} else {
if (max(sel) > d[4] || min(sel) < 1) {
note.AFNI(printf("sel outside of range 1:%d\n", d[4]))
return(NULL)
}
v <- brk[,,,sel]
d[4] <- length(sel)
dim(v) <- d
}
return(v)
}
dset.dimBRKarray <- function(dset) {
return(dimBRKarray(dset$brk))
}
index3Dto1D <- function(ii, dd) {
dd12 <- dd[1]*dd[2]
return( (ii[,1]-1)+
(ii[,2]-1)*dd[1]+
(ii[,3]-1)*dd12 + 1)
}
index1Dto3D <- function(ii, dd) {
dd12 <- dd[1]*dd[2]
ii <- ii-1
k <- ii %/% dd12
j <- ii %% dd12
i <- j %% dd[1]
j <- j %/% dd[1]
return(cbind(i+1,j+1,k+1))
}
dset.index3Dto1D <- function(ii, dset) {
dd <- dset$dim
return( index3Dto1D(ii,dd ))
}
dset.index1Dto3D <- function(ii, dset) {
dd <- dset$dim
return( index1Dto3D(ii,dd ))
}
dset.3DBRKarrayto1D <- function(dset) {
if (is.null(dset$brk)) {
err.AFNI("NULL or non existent dset$brk")
return(NULL);
}
if (dset.ndim(dset) == 1) return(dset)
dd <- dset.dimBRKarray(dset)
dim(dset$brk) <- c(prod(dd[1:3]),dd[4])
return(dset)
}
dset.1DBRKarrayto3D <- function(dset, dd=NULL) {
if (is.null(dset$brk)) {
err.AFNI("NULL or non existent dset$brk")
return(NULL);
}
if (dset.ndim(dset) == 3) return(dset)
ddi <- dset.dimBRKarray(dset)
if (is.null(dd)) { #get it from voxel header
dd <- dset$dim
}
if (prod(dd[1:3]) != prod(ddi[1])) {
err.AFNI("Dimension mismatch");
return(NULL)
}
dim(dset$brk) <- c(dd[1:3],ddi[2])
return(dset)
}
array2dset <- function (brk=NULL, format='BRIK', meth='AUTO') {
if (meth == 'AUTO') {
if (have_R_io()) meth <- 'clib' else meth <- 'Rlib'
}
if (is.vector(brk)) {
brk <- as.array(brk, dim=c(length(brk),1,1,1))
}
dd <- dimBRKarray(brk)
if (meth=='Rlib') {
z <- list(brk=brk,format=format, dim=dd,
delta=NULL, origin=NULL, orient=NULL)
class(z) <- "AFNI_R_dataset"
} else if (meth == 'clib') {
z <- list(brk=brk, NI_head<- list())
z$NI_head <- dset.attr(z$NI_head,
"IDCODE_STRING", val="SET_AT_WRITE_RANDOM")
z$NI_head <- dset.attr(z$NI_head, "ORIGIN", val=c(0,0,0), tp="float")
z$NI_head <- dset.attr(z$NI_head, "DELTA", val=c(1,1,1), tp="float")
z$NI_head <- dset.attr(z$NI_head, "ORIENT_SPECIFIC",
val = orcode.AFNI('RAI'))
if (!is.null(names(brk)))
z$NI_head <- dset.attr(z$NI_head, "BRICK_LABS",
val = names(brk));
class(z) <- "AFNI_c_dataset"
} else {
err.AFNI('bad meth');
return(NULL)
}
return(z)
}
typecode.AFNI <- function(typestr, def='MRI_short') {
if (is.null(typestr)) return(typecode.AFNI(def))
if (is.character(typestr)) {
if (typestr == 'MRI_byte' || typestr == 'byte') return(0)
if (typestr == 'MRI_short' || typestr == 'short') return(1)
if (typestr == 'MRI_int' || typestr == 'int') return(2)
if (typestr == 'MRI_float' || typestr == 'float') return(3)
if (typestr == 'MRI_double' || typestr == 'double') return(4)
if (typestr == 'MRI_complex' || typestr == 'complex') return(5)
if (typestr == 'MRI_rgb' || typestr == 'rgb') return(6)
err.AFNI(paste('Bad typecode ', typestr));
return(NULL);
} else {
if (typestr < 0 || typestr >6) {
err.AFNI(paste('Bad typecode ', typestr));
return(NULL);
}
return(typestr);
}
err.AFNI(paste('Should not be here ', typestr));
return(NULL);
}
orcode.AFNI <- function(orstr) {
if (is.character(orstr)) {
orcode <- c (-1,-1,-1)
orstr <- strsplit(orstr,'')[[1]]
for (i in 1:3) {
switch (tolower(orstr[i]),
r = orcode[i] <- 0,
l = orcode[i] <- 1,
p = orcode[i] <- 2,
a = orcode[i] <- 3,
i = orcode[i] <- 4,
s = orcode[i] <- 5)
if (orcode[i] < 0) {
err.AFNI(paste('Bad orientation code ', orstr));
return(NULL)
}
}
} else {
for (i in 1:3) {
orcode <- orstr
if (length(orcode) != 3 || orcode[i] < 0 || orcode[i] > 5) {
err.AFNI(paste('Bad orientation code ', orstr));
return(NULL)
}
}
}
return(orcode);
}
#------------------------------------------------------------------
# Extracted (and modified) from fmri library by Karsten Tabelow,
# tabelow@wias-berlin.de and Joerg Polzehl (polzehl@wias-berlin.de)
#
# Updates by ZSS & GC
#------------------------------------------------------------------
read.AFNI <- function(filename, verb = 0, ApplyScale = 1, PercMask=0.0,
meth = 'AUTO', forcedset = FALSE) {
if (meth == 'AUTO') {
if (have_R_io()) meth <- 'clib' else meth <- 'Rlib'
}
an <- parse.AFNI.name(filename);
if (verb > 1) {
show.AFNI.name(an);
}
if ( !forcedset && ((an$type == "1D" && an$ext != ".1D.dset")
|| an$type == "Rs" || an$type == "1Ds")) {
brk <- read.AFNI.matrix(an$orig_name)
z <- array2dset(brk, format=an$type)
return(z)
}
if (meth == 'clib' || an$type == 'NIML' ||
(an$type == "1D" && an$ext == ".1D.dset")) {
return(read.c.AFNI(filename, verb = verb, ApplyScale = 1, PercMask=0.0))
}
#If you have any selectors, use 3dbucket to get what you want, then read
#temp dset. This is an ugly fix for now, but will change it later if
#I/O is issue
if (verb) { cat ('Need tmp?\n'); }
rmtmp <- 0;
if (!is.null(an$brsel) || !is.null(an$rosel) || !is.null(an$rasel)) {
rmtmp <- 1;
#Changed below from >& /dev/null to > /dev/null 2>&1
#because system uses sh not tcsh
#Tx to G. Pagnoni & Co.
com <- paste ('3dcalc -overwrite -prefix ___R.read.AFNI.' ,
basename(an$pprefix),
' -a "', filename,'" -expr "a" > /dev/null 2>&1',
sep = '');
if (try(system(com)) != 0) {
warning(paste("Failed to execute:\n ", com),
immediate. = TRUE);
return(NULL);
}
an$pprefix <- paste('___R.read.AFNI.',basename(an$pprefix), sep = '');
if (!(exists.AFNI.name(head.AFNI.name(an)))) {
warning(paste("Failed to create: ",
head.AFNI.name(an), brik.AFNI.name(an), '\n'),
immediate. = TRUE);
return(NULL);
}
}
if (verb) { cat ('Checking existence\n'); }
if (!(exists.AFNI.name(head.AFNI.name(an)))) {
err.AFNI(paste("Failed to read: ", head.AFNI.name(an),
brik.AFNI.name(an)));
return(NULL);
}
#Cannot read compressed stuff (see size usage below)
if (verb) { cat ('Uncompressing\n'); }
uncompress.AFNI(head.AFNI.name(an));
conhead <- file(head.AFNI.name(an),"r")
header <- readLines(conhead)
close(conhead)
types <- NULL
args <- NULL
counts <- NULL
values <- NULL
for (i in 1:length(header)) {
if (regexpr("^type *= *", header[i]) != -1) {
tmptype <- strsplit(header[i]," *= *")[[1]][2]
types <- c(types,tmptype)
args <- c(args,strsplit(header[i+1]," *= *")[[1]][2])
tmpcounts <- as.numeric(strsplit(header[i+2]," *= *")[[1]][2])
counts <- c(counts,tmpcounts)
i <- i+3
tmpvalue <- ""
while ((regexpr("^$", header[i]) == -1) && (i <= length(header))) {
tmpvalue <- paste(tmpvalue,header[i])
i <- i+1
}
tmpvalue <- sub("^ +","",tmpvalue)
if ((tmptype == "integer-attribute") || (tmptype == "float-attribute")) {
tmpvalue <- as.numeric(strsplit(tmpvalue," +")[[1]])
}
values <- c(values,list(value=tmpvalue))
}
}
names(values) <- args
dx <- values$DATASET_DIMENSIONS[1]
dy <- values$DATASET_DIMENSIONS[2]
dz <- values$DATASET_DIMENSIONS[3]
dt <- values$DATASET_RANK[2]
scale <- values$BRICK_FLOAT_FACS
size <- file.info(brik.AFNI.name(an))$size/(dx*dy*dz*dt)
if (is.na(size)) {
err.AFNI("Failed to determine file size");
return(NULL);
}
if (regexpr("MSB",values$BYTEORDER_STRING[1]) != -1) {
endian <- "big"
} else {
endian <- "little"
}
if (min(abs(values$DELTA)) != 0) {
weights <-
abs(values$DELTA/min(abs(values$DELTA)))
} else {
weights <- NULL
}
# browser()
if (verb) { cat ('Reading Bin\n'); }
if (as.integer(size) == size) {
conbrik <- file(brik.AFNI.name(an),"rb")
# modified below by GC 12/2/2008
if (all(values$BRICK_TYPES==0) | all(values$BRICK_TYPES==1)) {
mybrk<- readBin( conbrik, "int", n=dx*dy*dz*dt, size=size,
signed=TRUE, endian=endian)
} else if (all(values$BRICK_TYPES==3)) {
mybrk<- readBin(conbrik, "numeric", n=dx*dy*dz*dt, size=size,
signed=TRUE, endian=endian) # float
} else {
err.AFNI("Cannot read datasets of multiple data types");
close(conbrik)
return(NULL);
}
close(conbrik)
dim(mybrk) <- c(dx,dy,dz,dt)
if (ApplyScale) {
if (verb) { cat ('Scaling\n'); }
#After this operation, size of mytt doubles if initially read as int
for (k in 1:dt) if (scale[k] != 0) mybrk[,,,k] <- scale[k] * mybrk[,,,k]
}
mask=NULL;
if (PercMask > 0.0) { #ZSS: Dunno what that is for.
# 0.75 was default for PercMask
mask <- array(TRUE,c(dx,dy,dz))
mask[mybrk[,,,1] < quantile(mybrk[,,,1],PercMask)] <- FALSE
}
z <- list(brk=mybrk,format=an$type,delta=values$DELTA,
origin=values$ORIGIN,
orient=values$ORIENT_SPECIFIC,
dim=c(dx,dy,dz,dt),weights=weights, header=values,mask=mask)
} else {
warning("Error reading file: Could not detect size per voxel\n")
z <- list(brk=NULL,format=an$type,delta=NULL,
origin=NULL,orient=NULL,dim=NULL,weights=NULL,
header=values,mask=NULL)
}
class(z) <- "AFNI_R_dataset"
attr(z,"file") <- paste(filename,"BRIK",sep="")
if (rmtmp == 1) {
if (verb) {
cat ('ZSS: Will remove tmp files\n');
}
#Changed below from >& to 2&>1 because system uses sh not tcsh
#Tx to G. Pagnoni & Co.
system('\\rm -f ___R.read.AFNI.* > /dev/null 2>&1');
}else{
if (verb) {
cat ('ZSS: No temps to remove\n');
}
}
invisible(z);
}
read.c.AFNI <- function(filename, verb = 0, ApplyScale = 1, PercMask=0.0) {
if (!is.loaded('R_THD_load_dset')) {
err.AFNI("Missing R_io.so");
return(NULL);
}
if (ApplyScale != 1 || PercMask != 0.0) {
err.AFNI(paste("This function is not ready to abide by ApplyScale != 1 or",
"PercMask != 0.0"));
return(NULL);
}
uopts = list( debug=verb );
rs2 <- .Call("R_THD_load_dset", name = as.character(filename), opts = uopts )
if (is.null(rs2)) return(NULL);
hatr <- parse.c.AFNI.head(rs2$head);
rs2$head <- NULL
ddb <- dset.attr(hatr, "dim")
dim(rs2$brk) <- ddb
rs2[['format']] <- 1
rs2[['delta']] <- dset.attr(hatr, "delta")
rs2[['origin']] <- dset.attr(hatr, "origin")
rs2[['orient']] <- dset.attr(hatr, "orient")
rs2[['dim']] <- dset.attr(hatr, "dim")
rs2[['NI_head']] <- hatr
rs2[['header']] <- NULL
rs2[['mask']] <- NULL
rs2[['weights']] <- NULL
class(rs2) <- "AFNI_c_dataset"
invisible(rs2)
}
parse.c.AFNI.head <- function (head) {
j <- 1
nel <- NULL;
nms <- NULL;
while (j <= length(head)) {
nelc <- r.NI_read_str_element(strsplit(head[j], "\n")[[1]], FALSE)
if (!is.null(nelc)) {
nn <- r.NI_get_attribute(nelc, "atr_name");
if (!is.null(nn)) {
nms <- c(nms, nn)
nel <- c(nel, list(nelc));
} else {
warn.AFNI("Have element that has no atr_name");
}
}
j <- j + 1;
}
if (!is.null(nel)) {
names(nel) <- nms
}
return(nel)
}
unparse.c.AFNI.head <- function (nel) {
j <- 1
head <- NULL;
nms <- NULL;
while (j <= length(nel)) {
#r.NI_write_str_element(nel[[j]],FALSE)
nm <- r.NI_get_attribute(nel[[j]],"atr_name")
if (!is.null(nm)) nms <- c(nms,nm)
else nms <- c(nms, 'noname');
head <- c(head, r.NI_write_str_element(nel[[j]],FALSE));
j <- j + 1;
}
if (length(nms)) names(head)<-nms
return(head)
}
#A funtion to create an AFNI header string
AFNIheaderpart <- function(type, name, value) {
a <- "\n"
a <- paste(a, "type = ", type, "\n", sep="")
a <- paste(a, "name = ", name, "\n", sep="")
if (regexpr("string",type) == 1) {
value <- paste("'", value, "~", sep="")
a <- paste(a, "count = ", nchar(value) - 1, "\n", sep ="")
a <- paste(a, value, "\n", sep="")
} else {
a <- paste(a, "count = ", length(value), "\n", sep ="")
j <- 0
while (j<length(value)) {
left <- length(value) - j
if (left>4) left <- 5
a <- paste(a, paste(value[(j+1):(j+left)],collapse=" "), "\n", sep=" ")
j <- j+5
}
}
return(a);
}
#Calculate the min, and max of BRIK data y
minmax <- function(y) {
r <- NULL;
for (k in 1:dim(y)[4]) {
r <- c(r,min(y[,,,k]),max(y[,,,k]))
};
return(r);
}
newid.AFNI.old <- function(ext=0) {
if (ext) { #in house
return(
paste('GCR_',paste(
sample(c(rep(0:9,each=5),LETTERS,letters),22,replace=TRUE),
collapse=''),
sep='') )
} else { #Call AFNI program
return(system('3dnewid -fun', intern=TRUE) );
}
}
AFNI.view2viewtype <- function(view="+orig"){
if (is.null(view)) return(0)
if (length(grep("orig", view, ignore.case=TRUE))) return(0)
if (length(grep("acpc", view, ignore.case=TRUE))) return(1)
if (length(grep("tlrc", view, ignore.case=TRUE))) return(2)
return(0)
}
AFNI.viewtype2view <- function(view=0){
if (is.null(view)) return("+orig")
if (view == 0) return("+orig")
if (view == 1) return("+acpc")
if (view == 2) return("+tlrc")
return("+orig")
}
dset.view <- function(dset) {
dd <- dset.attr(dset$NI_head, "SCENE_DATA")
return( AFNI.viewtype2view(dd[1]))
}
dset.gridmatch <- function(mdset, idset) {
if (is.null(mdset) || is.null(idset)) return(FALSE)
if (dset.nvox(mdset) != dset.nvox(idset)) return(FALSE)
return(TRUE);
}
write.c.AFNI <- function( filename, dset=NULL, label=NULL, space=NULL,
note=NULL, origin=NULL, delta=NULL,
orient=NULL,
idcode=NULL, defhead=NULL,
verb = 1,
maskinf=0, scale = TRUE,
overwrite=FALSE, addFDR=0,
statsym=NULL, view=NULL,
com_hist=NULL, TR=NULL, type=NULL) {
an <- parse.AFNI.name(filename);
if (verb > 1) {
show.AFNI.name(an);
}
if (class(dset) == "AFNI_c_dataset") {
#Assume all is in dset already
if (is.null(dset$NI_head)) {
err.AFNI("Old dset?");
return(NULL);
}
} else if (class(dset) == "AFNI_R_dataset") {
err.AFNI("Should not be with R dataset here")
return(NULL);
} else {
#Assume user sent in an array, not a dset
dset <- array2dset(dset,meth='clib')
}
if (!is.loaded('R_THD_write_dset')) {
err.AFNI("Missing R_io.so.");
return(NULL);
}
if (is.null(view)) {
if (is.null(defhead)) {
view <- "+tlrc"
}else {
view <- dset.view(defhead)
}
}
#Revert to old default if view is still not set
if (is.null(view)) view <- "+tlrc"
#Setup the basics (might be repetitive, oh well)
dset$NI_head <- dset.attr(dset$NI_head, "dim", val = dset.dimBRKarray(dset))
#Make sure TYPESTRING and SCENE_DATA are coherent
dset$NI_head <- dset.attr( dset$NI_head, "TYPESTRING",
default = '3DIM_GEN_FUNC')
dset$NI_head <- dset.attr( dset$NI_head, "SCENE_DATA",
default = c(AFNI.view2viewtype(view),11,3, rep(0,5)))
taxnums <- NULL;
taxfloats <- NULL;
taxoff <- NULL;
if (!is.null(defhead)) {
taxnums <- dset.attr( defhead, "TAXIS_NUMS");
taxfloats <- dset.attr( defhead, "TAXIS_FLOATS");
taxoff <- dset.attr( defhead, "TAXIS_OFFSETS");
}
#The user options
if (is.null(idcode)) idcode <- newid.AFNI();
if (!is.null(idcode)) dset$NI_head <-
dset.attr(dset$NI_head, "IDCODE_STRING", val=idcode)
if (!is.null(TR)) {
if (is.null(taxnums)) {
#77002 for seconds
taxnums <- c(dset.nval(dset), 0, 77002, -999, -999, -999, -999, -999);
}
if (is.null(taxfloats)) {
taxfloats <- c(0, TR, 0, 0, 0, -999999.000, -999999.000, -999999.0000);
}
if (is.null(taxoff)) {
taxoff <- rep(0, dset.nslize(dset));
}
dset$NI_head <- dset.attr( dset$NI_head, "TAXIS_NUMS", val = taxnums);
dset$NI_head <- dset.attr( dset$NI_head, "TAXIS_FLOATS", val = taxfloats);
dset$NI_head <- dset.attr( dset$NI_head, "TAXIS_OFFSETS", val = taxoff);
dset$NI_head <- dset.attr(dset$NI_head, "TR", val = TR);
} else { #use defhead values, setup from above
if (!is.null(taxnums)) {
dset$NI_head <- dset.attr( dset$NI_head, "TAXIS_NUMS", val = taxnums);
}
if (!is.null(taxfloats)) {
dset$NI_head <- dset.attr( dset$NI_head, "TAXIS_FLOATS",
val = taxfloats);
}
if (!is.null(taxoff)) {
dset$NI_head <- dset.attr( dset$NI_head, "TAXIS_OFFSETS",
val = taxoff);
}
}
if (is.null(origin) && !is.null(defhead))
origin <- dset.attr(defhead,"ORIGIN")
if (!is.null(origin)) dset$NI_head <-
dset.attr(dset$NI_head, "ORIGIN", val=origin)
if (is.null(delta) && !is.null(defhead))
delta <- dset.attr(defhead,"DELTA")
if (!is.null(delta)) dset$NI_head <-
dset.attr(dset$NI_head, "DELTA", val=delta)
IJK_TO_DICOM_REAL <- NULL
if (!is.null(defhead))
IJK_TO_DICOM_REAL <- dset.attr(defhead$NI_head,"IJK_TO_DICOM_REAL")
if (!is.null(IJK_TO_DICOM_REAL)) dset$NI_head <-
dset.attr(dset$NI_head, "IJK_TO_DICOM_REAL", val=IJK_TO_DICOM_REAL)
if (is.null(orient) && !is.null(defhead))
orient <- dset.attr(defhead,"ORIENT_SPECIFIC")
if (!is.null(orient)) dset$NI_head <-
dset.attr(dset$NI_head, "ORIENT_SPECIFIC",
val = orcode.AFNI(orient))
if (is.null(label) && !is.null(defhead))
label <- dset.attr(defhead,"BRICK_LABS")
if (!is.null(label)) dset$NI_head <-
dset.attr(dset$NI_head, "BRICK_LABS",
val = label);
if (is.null(space) && !is.null(defhead))
space <- dset.attr(defhead,"TEMPLATE_SPACE")
if (!is.null(space)) dset$NI_head <-
dset.attr(dset$NI_head, "TEMPLATE_SPACE",
val = space);
if (is.null(type) && !is.null(defhead)) {
tps = dset.attr(defhead,"BRICK_TYPES");
if (!is.null(tps) && length(tps)>0) type <- tps[1]
}
if (!is.null(type)) {
if (!is.null(type <- typecode.AFNI(type))) {
type <- rep(type, dset.dimBRKarray(dset)[4])
dset$NI_head <- dset.attr(dset$NI_head, "BRICK_TYPES", val=type)
}
}
if (is.null(note) && !is.null(defhead))
note <- dset.attr(defhead,"note")
if (!is.null(note)) dset$NI_head <-
dset.attr(dset$NI_head, "note",
val = note);
if (is.null(statsym) && !is.null(defhead))
statsym <- dset.attr(defhead,"statsym")
if (!is.null(statsym)) dset$NI_head <-
dset.attr(dset$NI_head, "statsym",
val = statsym);
if (maskinf) {
if (verb) note.AFNI("Masking infs");
dset$brk[!is.finite(dset$brk)]=0
}
uopts = list( scale = scale, overwrite=overwrite, debug=verb,
addFDR=addFDR, hist=com_hist);
dset$head <- unparse.c.AFNI.head(dset$NI_head)
rso <- .Call("R_THD_write_dset", name = as.character(filename),
dset = as.list(dset), opts = uopts)
invisible(dset);
}
write.AFNI <- function( filename, brk=NULL, label=NULL, space=NULL,
note=NULL, origin=NULL, delta=NULL,
orient=NULL,
idcode=NULL, defhead=NULL,
verb = 0,
maskinf=0, scale = TRUE,
meth='AUTO', addFDR=FALSE,
statsym=NULL, view=NULL,
com_hist=NULL, type=NULL, TR = NULL,
overwrite=FALSE) {
if (meth == 'AUTO') {
if (class(brk) == 'AFNI_R_dataset') {
meth = 'Rlib'
} else if (class(brk) == 'AFNI_c_dataset') {
meth = 'clib'
} else if (have_R_io()) {
meth <- 'clib'
} else {
meth <- 'Rlib'
}
}
an <- parse.AFNI.name(filename);
if (verb>1) {
show.AFNI.name(an);
}
if (an$type == "1D" || an$type == "Rs" || an$type == "1Ds") {
return(write.AFNI.matrix(drop(brk), an$orig_name))
}
if (meth == 'clib' || an$type == 'NIML' || class(brk) == "AFNI_c_dataset") {
return(write.c.AFNI(filename, dset=brk, label=label, space=space,
note=note, origin=origin, delta=delta,
orient=orient,
idcode=idcode, defhead=defhead,
verb = verb,
maskinf=maskinf, scale = scale, addFDR=addFDR,
statsym=statsym, view=view, com_hist=com_hist,
type=type, TR=TR, overwrite=overwrite))
}
#Switch to old default for view if not using write.c.AFNI
if (is.null(view)) view <- "+tlrc"
if (is.na(an$view)) {
err.AFNI('Bad filename for old writing method. Need view');
show.AFNI.name(an);
return(0)
}
if (addFDR && verb) {
warn.AFNI("addFDR is not supported for AFNI_R_dataset structs");
}
if (!is.null(statsym) && verb) {
warn.AFNI("statsym is not supported for AFNI_R_dataset structs");
}
#Call with brk if dset is sent in
if (class(brk) == "AFNI_R_dataset") {
return(write.AFNI( filename, brk$brk, label, space, note, origin, delta,
orient, idcode, defhead, verb, maskinf, scale))
}
#Set the defaults.
if (is.null(note)) note <- '';
if (is.null(idcode)) idcode <- newid.AFNI();
if (is.null(brk)) {
err.AFNI("NULL input");
return(0)
}
#make sure we don't have 3 dims array, need always 4d
ddb <- dimBRKarray(brk)
dim(brk) <- ddb
if (is.null(defhead)) { # No default header
if (verb) note.AFNI("Have no default header");
if (is.null(label)) label <- paste(c(1:ddb[4]),collapse='~');
if (is.null(space)) space <-'TLRC';
if (is.null(origin)) origin <- c(0,0,0)
if (is.null(delta)) delta <- c(4,4,4)
if (is.null(orient)) orient <- 'RAI'
} else { #When possible, call on default header
if (verb) note.AFNI("Have default header");
if (is.null(label)) {
if (!is.null(defhead$BRICK_LABS)) {
label <- gsub("^'", '', defhead$BRICK_LABS);
} else {
label <- paste(c(1:ddb[4]),collapse='~');
}
}
if (is.null(space)) {
if (!is.null(defhead$NI_head$TEMPLATE_SPACE$dat)) {
space <- defheadTEMPLATE_SPACE$dat;
} else {
space <- 'TLRC';
}
}
if (is.null(origin)) {
if (!is.null(defhead$ORIGIN)) {
origin <- defhead$ORIGIN;
} else {
origin <- c(0,0,0);
}
}
if (is.null(delta)) {
if (!is.null(defhead$DELTA)) {
delta <- defhead$DELTA;
} else {
delta <- c(4,4,4);
}
}
if (is.null(orient)) {
if (!is.null(defhead$ORIENT_SPECIFIC)) {
orient <- defhead$ORIENT_SPECIFIC;
} else {
orient <- 'RAI';
}
}
}
if (is.null(brk) || !is.numeric(brk)) {
err.AFNI("data array improperly formatted");
return(0);
}
if (verb) {
note.AFNI("Writing header")
}
conhead <- file(head.AFNI.name(an), "w")
writeChar(AFNIheaderpart("string-attribute","TEMPLATE_SPACE",space),
conhead,eos=NULL)
writeChar(AFNIheaderpart("string-attribute","HISTORY_NOTE",note),
conhead,eos=NULL)
writeChar(AFNIheaderpart("string-attribute","TYPESTRING","3DIM_HEAD_FUNC"),
conhead,eos=NULL)
writeChar(AFNIheaderpart("string-attribute","IDCODE_STRING",''),
conhead,eos=NULL)
writeChar(AFNIheaderpart("string-attribute","IDCODE_DATE",date()),
conhead,eos=NULL)
if (an$view == '+orig') { sv <- 0 }
else if (an$view == '+acpc') { sv <- 1 }
else if (an$view == '+tlrc') { sv <- 2 }
else { sv <- 0 }
writeChar(AFNIheaderpart("integer-attribute","SCENE_DATA",
c(sv,11,1,-999,-999,-999,-999,-999)),
conhead,eos=NULL)
writeChar(AFNIheaderpart("integer-attribute","ORIENT_SPECIFIC",
orcode.AFNI(orient)),
conhead,eos=NULL)
writeChar(AFNIheaderpart("float-attribute","ORIGIN",origin),
conhead,eos=NULL)
writeChar(AFNIheaderpart("float-attribute","DELTA",delta),
conhead,eos=NULL)
if (maskinf) {
if (verb) note.AFNI("Masking infs");
brk[!is.finite(brk)]=0
}
mm <- minmax(brk)
if (verb) {
note.AFNI("minmax");
str(mm)
}
writeChar(AFNIheaderpart("float-attribute","BRICK_STATS",mm),
conhead,eos=NULL)
writeChar(AFNIheaderpart("integer-attribute","DATASET_RANK",
c(3,ddb[4],0,0,0,0,0,0)),
conhead,eos=NULL)
writeChar(AFNIheaderpart("integer-attribute","DATASET_DIMENSIONS",
c(ddb[1:3],0,0)),
conhead,eos=NULL)
if (is.null(type)) type <- 1
writeChar(AFNIheaderpart("integer-attribute","BRICK_TYPES",
rep(typecode.AFNI(type),ddb[4])),
conhead,eos=NULL)
if (scale) {
if (verb) {
note.AFNI("Computing scaling factors");
}
scale_fac <- rep(0,ddb[4])
for (k in 1:ddb[4]) {
scale_fac[k] <- max(abs(mm[2*k-1]),abs(mm[2*k]))/32767
}
} else {
scale_fac <- rep(0,ddb[4])
}
if (verb) {
note.AFNI("scaling factors:");
str(scale_fac);
}
writeChar(AFNIheaderpart("float-attribute","BRICK_FLOAT_FACS",scale_fac),
conhead,eos=NULL)
writeChar(AFNIheaderpart("string-attribute","BRICK_LABS",
paste(label,collapse="~")),
conhead,eos=NULL)
writeChar(AFNIheaderpart("string-attribute","BYTEORDER_STRING","MSB_FIRST"),
conhead,eos=NULL)
#The tross logger
trlog <- Sys.getenv("AFNI_HISTDB_SCRIPT")
if ( trlog != '') {
writeChar(AFNIheaderpart("string-attribute","HISTDB_SCRIPT", trlog ),
conhead,eos=NULL) ;
}
close(conhead)
# Write BRIK
conbrik <- file(brik.AFNI.name(an), "wb")
if (0) { #ZSS: old method,
# runs out of memory for large dsets
# when scaling. Code kept here for testing
if (scale) {
if (verb) {
note.AFNI("Applying scaling ")
}
for (k in 1:ddb[4]) {
brk[,,,k] <- brk[,,,k] / scale_fac[k]
}
}
dim(brk) <- NULL #Don't know why this was done here ....
if (verb) {
note.AFNI("Writing brik")
}
writeBin(as.integer(brk), conbrik,size=2, endian="big")
} else {
if (verb) {
note.AFNI("Writing /Scaling, meth 2")
}
#Write on sub-brick at a time to reduce memory use.
# as.integer will allocate a new copy
if (scale) {
for (k in 1:ddb[4]) {
writeBin(as.integer(brk[,,,k] / scale_fac[k]),
conbrik,size=2, endian="big")
}
} else {
for (k in 1:ddb[4]) {
writeBin(as.integer(brk[,,,k]),
conbrik,size=2, endian="big")
}
}
}
close(conbrik)
return(1);
}
read.NIFTI <- function(filename) {
fileparts <- strsplit(filename,"\\.")[[1]]
ext <- tolower(fileparts[length(fileparts)])
if (ext == "nii") {
filename.nii <- filename
filename.hdr <- paste(c(fileparts[-length(fileparts)],"hdr"),collapse=".")
filename.img <- paste(c(fileparts[-length(fileparts)],"img"),collapse=".")
} else if (ext == "hdr") {
filename.hdr <- filename
filename.img <- paste(c(fileparts[-length(fileparts)],"img"),collapse=".")
} else if (ext == "img") {
filename.hdr <- paste(c(fileparts[-length(fileparts)],"hdr"),collapse=".")
filename.img <- filename
} else {
filename.nii <- paste(filename,".nii",sep="")
filename.hdr <- paste(filename,".hdr",sep="")
filename.img <- paste(filename,".img",sep="")
}
if ((ext != "hdr") && (ext != "img") && (!is.na(file.info(filename.nii)$size))) {
con <- file(filename.nii,"rb")
header <- read.NIFTI.header(con)
if (!(header$magic == "n+1") && !(header$magic == "ni1"))
warning("Hmmm! Dont see the magic NIFTI string! Try to proceed, but maybe some weird results will occur!");
bytes <- header$voxoffset - 348
header$extension <- readBin(con,"raw",bytes)
} else {
if (is.na(file.info(filename.hdr)$size) | (file.info(filename.hdr)$size < 348))
stop("Hmmm! This does not seem to be a NIFTI header (hdr/img-pair)! Wrong size or does not exist!");
con <- file(filename.hdr,"rb")
header <- read.NIFTI.header(con)
header$extension <- NULL
close(con)
if (is.na(file.info(filename.img)$size))
stop("Hmmm! This does not seem to be a NIFTI header (hdr/img-pair)! img-file not found!");
con <- file(filename.img,"rb")
}
dx <- header$dimension[2]
dy <- header$dimension[3]
dz <- header$dimension[4]
dt <- header$dimension[5]
endian <- header$endian
if (header$datatype == 1) { # logical
what <- "raw"
signed <- TRUE
size <- 1
} else if (header$datatype == 2) { # unsigned char????
what <- "int"
signed <- FALSE
size <- if (header$bitpix) header$bitpix/8 else 2
} else if (header$datatype == 4) { # signed short
what <- "int"
signed <- TRUE
size <- if (header$bitpix) header$bitpix/8 else 2
} else if (header$datatype == 8) { # signed integer
what <- "int"
signed <- TRUE
size <- if (header$bitpix) header$bitpix/8 else 4
} else if (header$datatype == 16) { # float
what <- "double"
signed <- TRUE
size <- if (header$bitpix) header$bitpix/8 else 4
} else if (header$datatype == 32) { # complex
what <- "complex"
signed <- TRUE
size <- if (header$bitpix) header$bitpix/8 else 8
} else if (header$datatype == 64) { # double
what <- "double"
signed <- TRUE
size <- if (header$bitpix) header$bitpix/8 else 8
} else { # all other
what <- "raw"
signed <- TRUE
size <- 1
}
brk <- readBin(con, what, n=dx*dy*dz*dt, size=size, signed=signed, endian=endian)
close(con)
if (min(abs(header$pixdim[2:4])) != 0) {
weights <-
abs(header$pixdim[2:4]/min(abs(header$pixdim[2:4])))
} else {
weights <- NULL
}
dim(brk) <- c(dx,dy,dz,dt)
mask <- array(TRUE,c(dx,dy,dz))
mask[brk[,,,1] < quantile(brk[,,,1],0.75)] <- FALSE
z <- list(brk=writeBin(as.numeric(brk),raw(),4),format="NIFTI",delta=header$pixdim[2:4],
origin=NULL,orient=NULL,dim=header$dimension[2:5],weights=weights,header=header,mask=mask)
class(z) <- "AFNI_R_dataset"
invisible(z)
}
extract.data <- function(z,what="data") {
if (!("AFNI_R_dataset"%in%class(z)) || !("AFNI_c_dataset"%in%class(z))) {
warning("extract.data: data not of class <AFNI_[cR]_dataset>. Try to proceed but strange things may happen")
}
if (what=="residuals") {
if(!is.null(z$resscale)){
brk <- readBin(z$res,"integer",prod(z$dim),2)*z$resscale
dim(brk) <- z$dim
} else {
warning("extract.data: No residuals available, returning NULL")
brk <- NULL
}
} else {
brk <- readBin(z$brk,"numeric",prod(z$dim),4)
dim(brk) <- z$dim
}
invisible(brk)
}
# GC 01/12/2015: Used in 3dMVM.T and 3dLME.R to construct GLT and GLF stuff
#gl_Constr2 <- function(n_gl, code, lop) { # n_gl: number of tests: lop$num_glt or lop$num_glf; code: lop$glfCode
# if (n_gl > 0) {
# outList <- vector('list', 3)
# outList[[1]] <- vector('list', n_gl)
# outList[[2]] <- vector('list', n_gl)
# outList[[3]] <- vector('list', n_gl)
# for (n in 1:n_gl) { # assuming each GLT has one slope involved and placed last
# # if(length(lop$QV)==0) outList[[1]][[n]] <- glfConstr(code[[n]], lop$dataStr) else {
# if((length(lop$QV)==0) & is.na(lop$vVars)) outList[[1]][[n]] <- glfConstr(code[[n]], lop$dataStr) else {
# if((length(lop$QV)>0) & any(lop$QV %in% code[[n]])) {
# QVpos <- which(code[[n]] %in% lop$QV)
# if(is.na(code[[n]][QVpos+2])) { # test for covariate effect
# if(QVpos==1) outList[[1]][[n]] <- NA else
# outList[[1]][[n]] <-glfConstr(code[[n]][-c(QVpos, QVpos+1)], lop$dataStr)
# outList[[2]][[n]] <- code[[n]][QVpos]
# } else { # effect at a specific covariate value
# if(QVpos==1) outList[[1]][[n]] <- NA else
# outList[[1]][[n]] <-glfConstr(code[[n]][-(QVpos:(QVpos+2))], lop$dataStr)
# outList[[3]][[n]] <- as.numeric(code[[n]][QVpos+2])
# names(outList[[3]][[n]]) <- code[[n]][QVpos]
# } # if(is.na(lop$gltCode[[n]][QVpos+2]))
# } else if(!is.na(lop$vVars) & any(lop$vQV %in% code[[n]])) { # voxel-wise covariate
# vQVpos <- which(code[[n]] %in% lop$vQV)
# if(is.na(code[[n]][vQVpos+2])) { # test for covariate effect
# if(vQVpos==1) outList[[1]][[n]] <- NA else
# outList[[1]][[n]] <-glfConstr(code[[n]][-c(vQVpos, vQVpos+1)], lop$dataStr)
# outList[[2]][[n]] <- code[[n]][vQVpos]
# } else { # effect at a specific covariate value
# if(vQVpos==1) outList[[1]][[n]] <- NA else
# outList[[1]][[n]] <-glfConstr(code[[n]][-(vQVpos:(vQVpos+2))], lop$dataStr)
# outList[[3]][[n]] <- as.numeric(code[[n]][vQVpos+2])
# names(outList[[3]][[n]]) <- code[[n]][vQVpos]
# } # if(is.na(lop$gltCode[[n]][vQVpos+2]))
# } else outList[[1]][[n]] <- glfConstr(code[[n]], lop$dataStr) # if((length(lop$QV)>0) & any(lop$QV %in% lop$gltCode[[n]]))
# }
# }
# }
# return(outList)
#}
gl_Constr <- function(n_gl, code, lop) { # n_gl: number of tests: lop$num_glt or lop$num_glf; code: lop$glfCode
if (n_gl > 0) {
outList <- vector('list', 3)
outList[[1]] <- vector('list', n_gl)
outList[[2]] <- vector('list', n_gl)
outList[[3]] <- vector('list', n_gl)
comQV <- c(lop$QV, lop$vVars)
for (n in 1:n_gl) { # assuming each GLT has one slope involved and placed last
if(length(comQV)==0) outList[[1]][[n]] <- glfConstr(code[[n]], lop$dataStr) else {
if((length(comQV)>0) & any(comQV %in% code[[n]])) {
QVpos <- which(code[[n]] %in% comQV)
if(is.na(code[[n]][QVpos+2])) { # test for covariate effect
if(QVpos==1) outList[[1]][[n]] <- NA else
outList[[1]][[n]] <-glfConstr(code[[n]][-c(QVpos, QVpos+1)], lop$dataStr)
outList[[2]][[n]] <- code[[n]][QVpos]
} else { # effect at a specific covariate value
if(QVpos==1) outList[[1]][[n]] <- NA else
outList[[1]][[n]] <-glfConstr(code[[n]][-(QVpos:(QVpos+2))], lop$dataStr)
outList[[3]][[n]] <- as.numeric(code[[n]][QVpos+2])
names(outList[[3]][[n]]) <- code[[n]][QVpos]
} # if(is.na(lop$gltCode[[n]][QVpos+2]))
} else outList[[1]][[n]] <- glfConstr(code[[n]], lop$dataStr) # if((length(lop$QV)>0) & any(lop$QV %in% lop$gltCode[[n]]))
}
if(is.null(outList[[1]][[n]])) errex.AFNI(paste("Inappropriate coding in test No.", n, "! \n ", sep = ""))
}
}
return(outList)
}
beauchamplab/rave documentation built on June 18, 2024, 10:26 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
#------------------------------------------------------------------
# Global Variables
#------------------------------------------------------------------
BATCH_MODE <<- 0 #Initialize batch mode flag to 0
R_io <<- -1
SHOW_TRC <<- FALSE
#------------------------------------------------------------------
# Functions for library loading
#------------------------------------------------------------------
find.in.path <- function(file) { #Pretty much same as first.in.path
sys_path = strsplit(Sys.getenv('PATH'),':')[[1]]
ff <- paste(sys_path,'/', file, sep='')
ff = ff[file.exists(ff)]
ff = base::normalizePath(ff)
aa <- gsub('//','/',ff[1], fixed=TRUE)
if (is.na(aa)) aa <- NULL
return(aa)
}
#Return cat's file option to format the help output
#for TXT or sphinx purposes
help.cat.file.AFNI <- function (pname=NULL, targ='TXT') {
if (is.null(pname)) {
err.AFNI("NULL name for help function")
return("")
}
if (targ == 'TXT') {
dopt = '-hdoc_2_txt'
} else if (targ == 'SPX') {
dopt = '-hdoc_2_spx'
} else if (targ == 'ASPX') {
dopt = '-hdoc_2_aspx'
} else if (targ == 'RAW') {
return("")
} else {
warn.AFNI(paste("targ", targ,"unknown. Assuming 'TXT'"))
dopt = '-hdoc_2_txt'
}
return(paste("|apsearch ", dopt ,pname,"-"))
}
#print warnings a la AFNI
prompt.AFNI <- function (str='I await', choices=c('y','n'), vals=NULL) {
if (!is.null(vals) && length(vals) != length(choices)) {
err.AFNI(paste("Have ", length(choices), "options, but",
length(vals), "values to return"))
return(0)
}
choices[1]<-toupper(choices[1])
kk<-vector(length=0)
spr <- paste(str," [",paste(choices, collapse='|'),"]:", sep='')
if (BATCH_MODE) { #Only two choices available for this one!
if (length(choices)!=2) {
err.AFNI("This one can't run in batch mode for more than 2 choices")
return(0)
}
ss <- sys.AFNI(paste(
"prompt_user -pause '", spr,"'", sep='' ))
if (ss$out == "0") return(2) #The not default
else return(1) #The default
} else {
while (length(kk) == 0) {
cat(spr)
bb <- readLines(n=1)
if (bb == '') {
kk <- 1
} else {
kk <-which(tolower(choices) == tolower(bb))
}
}
if (!is.null(vals)) {
return(vals[kk])
} else {
return(kk)
}
}
return(0)
}
set.AFNI.msg.trace <- function (vv=FALSE) { SHOW_TRC <<- vv }
warn.AFNI <- function (str='Consider yourself warned',
callstr=NULL,
newline=TRUE) {
if (is.null(callstr)) {
if (SHOW_TRC) callstr <- who.called.me(TRUE)
else callstr <- ''
}
nnn<-''
if (newline) nnn <- '\n'
if (BATCH_MODE) ff <- stderr()
else ff <- ''
cat( '\n', 'oo Warning: ', callstr,'\n ',
paste(str, collapse=''), nnn,
sep='', file = ff)
}
err.AFNI <- function (str='Danger Danger Will Robinson',
callstr=NULL,
newline=TRUE) {
if (is.null(callstr)) {
if (SHOW_TRC) callstr <- who.called.me(TRUE)
else callstr <- ''
}
nnn<-''
if (newline) nnn <- '\n'
if (BATCH_MODE) ff <- stderr()
else ff <- ''
cat( '\n', '** Error: ', callstr,'\n ',
paste(str, collapse=''), nnn,
sep='', file = ff)
}
errex.AFNI <- function (str='Alas this must end',
callstr=NULL, newline=TRUE) {
if (is.null(callstr)) {
if (SHOW_TRC) callstr <- who.called.me(TRUE)
else callstr <- ''
}
err.AFNI(str,callstr, newline)
exit.AFNI(str='\n Execution halted',stat=1)
}
note.AFNI <- function (str='May I speak frankly?',
callstr=NULL, newline=TRUE, tic=1,
trimtrace=30) {
if (is.null(callstr)) {
if (SHOW_TRC) callstr <- who.called.me(TRUE, trim=trimtrace)
else callstr <- ''
}
nnn<-''
if (newline) nnn <- '\n'
if (BATCH_MODE) ff <- stderr()
else ff <- ''
if (tic == 1) {
tm <- format(Sys.time(), " @ %H:%M:%S")
} else if (tic==2) {
tm <- format(Sys.time(), " @ %a %b %d %H:%M:%S %Y")
} else tm <- ''
cat( '\n', '++ Note: ', callstr,
sprintf('%s\n ', tm),
paste(str, collapse=''),nnn,
sep='', file = ff)
}
exit.AFNI <- function(str='The piano has been drinking.', stat=0) {
if (BATCH_MODE) {
quit(save='no', status = stat)
} else {
#note.AFNI(str)
stop(str)
}
}
#Locate and load R_io.so
set_R_io <- function() {
rio <- 0
ll <- find.in.path('R_io.so')
if (!is.null(ll)) {
dd <- try(dyn.load(ll), silent=TRUE)
if (dd[[1]]!="R_io") {
warn.AFNI(paste("Failed to load R_io.so with this error message:\n"))
dyn.load(ll)
} else {
rio <- 1
}
}
return(rio)
}
if (R_io == -1) {
R_io <<- set_R_io()
}
have_R_io <- function() {
if (R_io == 1) return(TRUE) else return(FALSE)
}
libLoad <- function(myLib) {
sucLoad <- FALSE
sucCheck <- FALSE
try(sucLoad <- library(myLib, character.only = TRUE, logical.return = TRUE))
if (sucLoad) {
print(sprintf("Package %s successfully loaded!", myLib)); sucCheck <- TRUE
} else {
if (BATCH_MODE == 1) {
err.AFNI(paste( "Need to install package ",myLib,
"\n Start an interactive R session then run:\n",
"\n update.packages(checkBuilt=TRUE, ask=FALSE) ",
"#Not mandatory, but recommended",
"\n install.packages('",myLib,"')\n",
"\n Then quit R and rerun your command.\n",
sep=''))
} else {
try(install.packages(myLib))
try(sucLoad <- library(myLib, character.only = TRUE,
logical.return = TRUE))
if (sucLoad) print(sprintf("Package %s successfully loaded...", myLib))
}
}
return(sucLoad)
}
pkgLoad <- function(pkgs) {
# install packages not already loaded:
pkgs_miss <- pkgs[which(!pkgs %in% installed.packages()[, 1])]
if (length(pkgs_miss) > 0) {
install.packages(pkgs_miss, dep=TRUE, repos='http://watson.nci.nih.gov/cran_mirror/')
}
# load packages not already loaded:
attached <- search()
attached_pkgs <- attached[grepl("package", attached)]
need_to_attach <- pkgs[which(!pkgs %in% gsub("package:", "", attached_pkgs))]
if (length(need_to_attach) > 0) {
for (i in 1:length(need_to_attach)) if(require(need_to_attach[i], character.only = TRUE))
cat('Package ', need_to_attach[i], ' loaded successfully!\n\n', sep='')
}
if (length(need_to_attach) == 0) {
message("\n ...All required packages were already loaded!\n")
}
}
#------------------------------------------------------------------
# Functions to deal with AFNI file names
#------------------------------------------------------------------
strip.extension <- function (filename, extvec=NULL, verb=0) {
n <- list()
if (is.null(extvec)) {
ff <- strsplit(filename, '\\.')[[1]]
if (length(ff) > 1) {
n$ext <- paste('.',ff[length(ff)], sep='')
n$name_noext <- paste(ff[1:length(ff)-1],collapse='.')
} else {
n$ext <- ''
n$name_noext <- filename
}
} else {
n$ext <- ''
n$name_noext <- filename
for (ex in extvec) {
patt <- paste('\\',ex,'$',collapse='', sep='')
if (length(grep(patt, filename))) {
n$ext <- ex
n$name_noext <- sub(patt,'',filename)
return(n)
}
}
}
return(n)
}
parse.name <- function (filename, extvec=NULL, verb=0) {
n <- list()
ff <- strsplit(filename, .Platform$file.sep)[[1]]
n$name <- ff[length(ff)]
n$path <- '.'
if (length(ff)>1)
n$path <- paste(ff[1:length(ff)-1],collapse=.Platform$file.sep)
n$path <- paste(n$path, .Platform$file.sep, sep="")
n2 <- strip.extension(n$name, extvec, verb)
n$ext <- n2$ext
n$name_noext <- n2$name_noext
return(n)
}
is.AFNI.1D.string <- function(t) {
if (length(grep('^1D:',t))) return(TRUE)
return(FALSE)
}
eval.AFNI.string.help <- function() {
return("
Data Strings:
-------------
You can specify input matrices and vectors in a variety of
ways. The simplest is by specifying a .1D file with all
the trimmings of column and row selectors. You can also
specify a string that gets evaluated on the fly.
For example: '1D: 1 4 8' evaluates to a vector of values 1 4 and 8.
Also, you can use R expressions such as: 'R: seq(0,10,3)'
")
}
eval.AFNI.1D.string <- function (t, verb=0, nmax=0) {
#remove 1D:
t <- sub("^1D:","",t)
#any transpose?
doTr = FALSE
if (length(grep("'$",t))) {
t<-sub("'$",'',t)
doTr = TRUE
}
#replace commas with space
t<-gsub(",",' ',t)
#remove multiple blanks
t <- deblank.string(t, middle=TRUE)
vvf <- vector(length = 0, mode="numeric")
#replace .. with : and split at 'space'
s <- strsplit(sub("..",":",t, fixed=TRUE), " ")[[1]]
#replace $ with nmax if possible
s <- gsub("[$]",as.character(nmax),s)
#Now loop and form vector of components
for (si in s) {
if (verb) cat ("working ", si, "\n")
if (length(grep(":",si))) {
sss <- strsplit(si,'[(*)]')[[1]]
if (length(sss)>1) {
se <- sss[2]
si <- strsplit(sss[1],":")[[1]]
vv <- eval(parse(text=sprintf('seq(from=%s,to=%s,by=%s)',
si[1], si[2], se) ))
} else {
vv <- eval(parse(text=si))
}
} else if (length(grep("@",si))) {
ssi = as.numeric(strsplit(si,"@")[[1]])
#cat ("ssi = ", ssi, "\n")
vv = rep(ssi[2], ssi[1])
} else {
vv = as.numeric(si)
}
#cat(si," = ",vv, "\n")
vvf <- c(vvf, vv)
#cat("vvnow = ",vvf, "\n")
}
if (doTr) return(t(vvf))
else return(vvf)
}
is.AFNI.R.string <- function(t) {
if (length(grep('^R:',t))) return(TRUE)
return(FALSE)
}
eval.AFNI.R.string <- function (t) {
t <- sub("^R:","",t)
return(eval(parse(text=t)))
}
eval.AFNI.string <- function (t) {
if (is.AFNI.1D.string(t)) return(eval.AFNI.1D.string(t))
else if (is.AFNI.R.string(t)) return(eval.AFNI.R.string(t))
else return(NULL)
}
date.stamp <- function (fancy=FALSE) {
if (fancy) {
return(gsub(' ','_',date()))
} else {
return(format(Sys.time(), "%d%m%y-%H%M%S"))
}
}
parse.AFNI.name.selectors <- function(filename,verb=0) {
n <- list()
n$brsel<- NULL
n$rosel<- NULL
n$rasel<- NULL
n$insel<- NULL
selecs <- strsplit(filename,"\\[|\\{|<|#")[[1]]
n$name <- selecs[1]
for (ss in selecs[2:length(selecs)]) {
if (length(grep("]",ss))) {
n$brsel <- strsplit(ss,"\\]")[[1]][1]
} else if (length(grep("}",ss))) {
n$rosel <- strsplit(ss,"\\}")[[1]][1]
} else if (length(grep(">",ss))) {
n$rasel <- strsplit(ss,">")[[1]][1]
}
}
selecs <- strsplit(filename,"#")[[1]]
if (length(selecs) > 1) {
n$insel <- selecs[2]
}
return(n)
}
parse.AFNI.name <- function(filename, verb = 0) {
if (filename == '-self_test') { #Secret testing flag
note.AFNI('Function running in test mode')
show.AFNI.name(parse.AFNI.name('DePath/hello.DePrefix', verb))
show.AFNI.name(parse.AFNI.name('DePath/DePrefix+acpc', verb))
show.AFNI.name(parse.AFNI.name('DePath/DePrefix+acpc.', verb))
show.AFNI.name(parse.AFNI.name('DePath/DePrefix+acpc.HEAD', verb))
show.AFNI.name(parse.AFNI.name('DePath/DePrefix+acpc.BRIK.gz', verb))
show.AFNI.name(parse.AFNI.name('DePath/DePrefix+acpc.HEAD[23]', verb))
show.AFNI.name(
parse.AFNI.name('DePath/DePrefix+acpc.HEAD[DeLabel]{DeRow}', verb))
show.AFNI.name(
parse.AFNI.name('DePath/DePrefix+acpc[DeLabel]{DeRow}', verb))
show.AFNI.name(
parse.AFNI.name('DePath/DePrefix+acpc.[DeLabel]{DeRow}', verb))
return(NULL)
}
an <- list()
an$view <- NULL
an$pprefix <- NULL
an$brsel <- NULL
an$rosel <- NULL
an$rasel <- NULL
an$insel <- NULL
an$type <- NULL
an$path <- NULL
an$orig_name <- filename
an$file <- NULL
if (verb) { cat ('Parsing >>',filename,'<<\n', sep=''); }
if (!is.character(filename)) {
warning(paste('filename >>',
filename, '<< not a character string\n', sep=''),
immediate. = TRUE)
traceback()
return(NULL)
}
#Deal with special names:
if (length(grep("^1D:.*$",filename))) {
an$type = '1Ds'
return(an)
} else if (length(grep("^R:.*$",filename))) {
an$type = 'Rs'
return(an)
}
#Deal with selectors
n <- parse.AFNI.name.selectors(filename, verb)
filename <- n$name
an$file <- n$name
an$brsel <- n$brsel
an$rosel <- n$rosel
an$rasel <- n$rasel
an$insel <- n$insel
#Remove last dot if there
filename <- sub('\\.$','',filename)
#NIFTI?
n <- strip.extension(filename, c('.nii', '.nii.gz'), verb)
if (n$ext != '') {
an$ext <- n$ext
an$type <- 'NIFTI'
an$pprefix <- n$name_noext
} else {
#remove other extensions
n <- strip.extension(filename, c('.HEAD','.BRIK','.BRIK.gz',
'.BRIK.bz2','.BRIK.Z',
'.1D', '.1D.dset',
'.niml.dset',
'.' ),
verb)
if (n$ext == '.1D' || n$ext == '.1D.dset') {
an$type <- '1D'
} else if (n$ext == '.niml.dset') {
an$type <- 'NIML'
} else {
an$type <- 'BRIK'
}
if (n$ext == '.') {
n$ext <- ''
}
an$ext <- n$ext
filename <- n$name_noext
n <- strip.extension(filename, c('+orig','+tlrc','+acpc'), verb)
if (n$ext != '') {
an$view <- n$ext
} else {
an$view <- NA
}
an$pprefix <- n$name_noext
}
#a prefix with no path
an$prefix <- basename(an$pprefix)
#and the path
an$path <- dirname(an$orig_name)
if (verb > 2) {
note.AFNI("Browser not active")
# browser()
}
if ( an$type != '1D' && (
!is.null(an$brsel) || !is.null(an$rosel) ||
!is.null(an$rasel) || !is.null(an$insel))) {
#Remove trailing quote if any
an$prefix <- gsub("'$", '', an$prefix)
an$prefix <- gsub('"$', '', an$prefix)
an$pprefix <- gsub("'$",'', an$pprefix)
an$pprefix <- gsub('"$','', an$pprefix)
}
if ( an$type != 'BRIK' ) {
#Put the extension back on
an$pprefix <- paste(an$pprefix,an$ext, sep='')
an$prefix <- paste(an$prefix,an$ext, sep='')
}
return(an)
}
exists.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an)
ans <- 0
if (file.exists(head.AFNI.name(an))) ans <- ans + 1;
if (file.exists(brik.AFNI.name(an)) ||
file.exists(paste(brik.AFNI.name(an),'.gz', sep='')) ||
file.exists(paste(brik.AFNI.name(an),'.Z', sep=''))) ans <- ans + 2
return(ans)
}
#Read parameters for -dataTable options (e.g. 3dLME)
#If opts is just one string starting with '@' such as
#@PPP, the parameters are read from text file PPP and
#returned as if they were found on the command line
dataTable.AFNI.parse <- function(opts) {
if (is.null(opts) || length(opts) == 0) {
return(NULL)
}
if (length(opts) == 1 && length(grep('^@.*$',opts))) {
ff <- strsplit(opts, '')[[1]]
ff <- paste(ff[2:length(ff)], sep='', collapse='')
if (!file.exists(ff)) {
err.AFNI(paste("table file", ff, "does not exist"));
return(NULL)
}
#Can't read as table because users might have EOL \ from shell command
#opts <- scan(ff, what='character')
opts <- scan(ff, what='character', na.strings="") # GC 08/15/2014: modified to allow NA in dataTable
opts <- opts[grep('[^\\\\]',opts)]
}
return(opts);
}
used.AFNI.prefix <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
ans <- 0
ans <- exists.AFNI.name(an)
if (ans > 0) return(1)
if (exists.AFNI.name(sprintf('%s+orig', an$pprefix))) return(1)
if (exists.AFNI.name(sprintf('%s+tlrc', an$pprefix))) return(1)
if (exists.AFNI.name(sprintf('%s+acpc', an$pprefix))) return(1)
return(0)
}
pprefix.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
return(an$pprefix);
}
prefix.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
return(an$prefix);
}
view.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
return(an$view);
}
pv.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
return(paste(an$pprefix,an$view,sep=''));
}
head.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
if (an$type == 'BRIK' && !is.na(an$view)) {
return(paste(an$pprefix,an$view,".HEAD",sep=''));
} else {
return((an$orig_name));
}
}
brik.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
if (an$type == 'BRIK' && !is.na(an$view)) {
return(paste(an$pprefix,an$view,".BRIK",sep=''));
} else {
return((an$orig_name));
}
}
compressed.AFNI.name <- function(an) {
if (is.character(an)) an <- parse.AFNI.name(an);
if (length(grep('\\.gz$', an$ext))) {
return('gz')
} else if (length(grep('\\.bz2$', an$ext))) {
return('bz2')
} else if (length(grep('\\.Z$', an$ext))) {
return('Z')
} else {
return('')
}
}
modify.AFNI.name <- function (name, what="append", val="_new", cwd=NULL) {
if (!is.loaded('R_SUMA_ParseModifyName')) {
err.AFNI("Missing R_io.so");
return(NULL);
}
an <- .Call("R_SUMA_ParseModifyName",
name = name,
what = what,
val = val,
cwd = cwd)
return(an)
}
AFNI.command.history <- function(ExecName=NULL, args=NULL, ohist=NULL) {
if (!is.loaded('R_SUMA_HistString')) {
return(NULL);
}
an <- .Call("R_SUMA_HistString", ExecName, args, ohist)
return(an)
}
uncompress.AFNI <- function(an, verb = 1) {
if (is.character(an)) an <- parse.AFNI.name(an);
ans <- 0
zz <- paste(brik.AFNI.name(an),'.Z', sep='');
if (file.exists(zz)) {
if (verb) {
cat ('Uncompressing', zz, '\n');
}
system(paste('uncompress', zz));
}
zz <- paste(brik.AFNI.name(an),'.gz', sep='');
if (file.exists(zz)) {
if (verb) {
cat ('gzip -d-ing', zz, '\n');
}
system(paste('gzip -d', zz));
}
zz <- paste(brik.AFNI.name(an),'.bz2', sep='');
if (file.exists(zz)) {
if (verb) {
cat ('bzip2 -d-ing', zz, '\n');
}
system(paste('bzip2 -d', zz));
}
return(an);
}
copy.AFNI.dset <- function(an, ancp, overwrite=FALSE, verb = 1) {
if (is.character(an)) an <- parse.AFNI.name(an);
if (is.character(ancp)) ancp <- parse.AFNI.name(ancp);
ans <- 0
if (!exists.AFNI.name(an)) {
err.AFNI(sprintf('Dset %s not found', an$file));
return(ans)
}
if (exists.AFNI.name(ancp) && !overwrite) {
err.AFNI(sprintf('New dset %s already exists', ancp$file));
return(ans)
}
opt <- ''
if (overwrite) opt <- sprintf('%s %s', opt, -overwrite);
com = sprintf('3dcopy %s %s', an$file, ancp$file);
system(com, ignore.stderr = TRUE, intern=TRUE)
if (!exists.AFNI.name(ancp)) {
err.AFNI(sprintf('New dset %s not created', ancp$file));
return(ans)
}
return(1)
}
show.AFNI.name <- function(an) {
cat ('\n',
'uname=', an$orig_name, '\n',
'file=', an$file,'\n',
'type =', an$type,'\n',
'pref.=', an$pprefix, '\n',
'view =', an$view, '\n',
'head =', head.AFNI.name(an), '\n',
'brik =', brik.AFNI.name(an), '\n',
'brsel=', an$brsel, '\n',
'rosel=', an$rosel, '\n',
'rasel=', an$rasel, '\n',
'insel=', an$insel, '\n',
'compr=', compressed.AFNI.name(an), '\n',
'exist=', exists.AFNI.name(an), '\n');
}
#------------------------------------------------------------------
# Functions to parse command-line arguments
#------------------------------------------------------------------
value.AFNI.args <- function(name, ops) {
ifnd <- which(name == names(ops));
if (length(ifnd)) {
vv <- vector(typeof(ops[ifnd][[1]]));
for (i in 1:length(ifnd)) {
vv <- c (vv,ops[i][[1]])
}
return(vv);
}
return(NULL);
}
show.AFNI.args <- function (ops, verb=0, adieu=FALSE, hstr='') {
if (is.null(ops)) {
cat ('NULL options\n');
} else {
cat (hstr,'Allowed Options:\n');
if (length(ops[['allowed_options']])) {
for (i in 1:length(ops[['allowed_options']])) {
cat (' ', ops[['allowed_options']][i], '\n');
}
} else {
cat ('whatever grinds your beans');
}
if (verb) {
cat (hstr, 'User options:\n');
for (i in 1:length(ops)) {
if ((names(ops)[i] != 'allowed_options')) {
cat (' ', names(ops)[i], ': ', ops[[i]],'\n');
}
}
}
}
if (adieu) exit.AFNI(0)
}
check.AFNI.args <- function ( ops, params = NULL, verb=0) {
if (!is.null(params) && !is.null(ops)) {
for (i in 1:length(ops)) {
if (verb) {
str(names(ops)[i])
str(params[names(ops)[i]][[1]])
cat('\nChecking on ', paste(ops[[i]],collapse=','),'\n');
}
ipar <- which(names(ops)[i] == names(params));
if (length(ipar)) {
pp <- params[ipar[1]][[1]]['count'][[1]];
opsvec <- ops[[i]];
if (length(pp) == 1) { #exact number
if (length(opsvec) != pp) {
#browser()
msg <- paste( 'Expecting ',pp, ' parameters for option "',
names(ops)[i], '".\n Have ',
length(opsvec), ' parameter(s) in string "',
paste(opsvec, collapse = ' '),
'" instead.', sep = '')
if (length(opsvec) > pp &&
length(grep('^-[a-z,A-Z]', opsvec[1+pp]))) {
msg <- paste( msg, '\n NOTE that ', opsvec[1+pp],
' in bad option above is not a recognized option.',
collapse = '', sep = '' );
}
err.AFNI(msg);
return(0);
}
} else if (length(pp) == 2) { #range
if (length(opsvec) < pp[1] || length(opsvec) > pp[2]) {
if (pp[2] == Inf) {
msg <- paste( 'Expecting more than ',pp[1],
' parameters for option "',
names(ops)[i], '".\n Have ',
length(opsvec), ' parameter(s) in string "',
paste(opsvec, collapse = ' '),
'" instead.', sep = '')
if (length(opsvec) > 0 &&
length(grep('^-[a-z,A-Z]', opsvec[1]))) {
msg <- paste( msg, '\n NOTE that ', opsvec[1],
' in bad option above is not a recognized option.',
collapse = '', sep = '' );
}
err.AFNI(msg);
} else {
msg <- paste( 'Expecting ',pp[1], ' to ', pp[2],
' parameters for option "',
names(ops)[i], '".\n Have ',
length(opsvec), ' parameter(s) in string "',
paste(opsvec, collapse = ' '),
'" instead.', sep = '');
if (length(opsvec) > pp[2] &&
length(grep('^-[a-z,A-Z]', opsvec[1+pp[2]]))) {
msg <- paste( msg, '\n NOTE that ', opsvec[1+pp[2]],
' in bad option above is not a recognized option.',
collapse = '', sep = '' );
}
err.AFNI(msg);
}
return(0);
}
} else {
warning(paste( 'I do not know what to do here'),
immediate. = TRUE);
return(0);
}
} else {
#ok
}
}
}
return(1); #all ok
}
apl <- function ( n = 0, d=NA, h=NULL, dup=FALSE ) {
return(list('count'=n, 'default'=d, help=h, duplicate_ok=dup));
}
load.debug.AFNI.args <- function ( fnm = NULL) {
if (is.null(fnm)) {
fnm <- system('ls *.dbg.AFNI.args .*.dbg.AFNI.args', ignore.stderr = TRUE, intern=TRUE)
}
if (length(fnm) > 1) {
err.AFNI(paste("More than one argument file found:",
paste(fnm, collapse=' ', sep=' '), ".\n",
"Try one of:\n",
paste("load.debug.AFNI.args('",fnm,"')",
collapse='\n', sep=''),
sep = ''));
return(FALSE);
}else if (length(fnm) == 0) {
err.AFNI("No files to load");
return(FASLE);
}
load(fnm)
ok_auto <- FALSE
if (exists('rfile')) {
ss <- sprintf("\n Start interactive code with:\n\nsource('%s')\n",
rfile);
ok_auto <- TRUE
} else {
ss <-
sprintf("\n Start interactive code with:\n\nsource('yourRfile')\n");
}
if (exists('args')) {
note.AFNI(sprintf("Setting .DBG_args from args in %s%s", fnm, ss));
.DBG_args <<- args
} else {
err.AFNI(sprintf("Variable args not in %s", fnm));
return(FALSE)
}
if (ok_auto) {
if (prompt.AFNI("Execute source command?", c('y','n')) == 1) {
source(rfile)
}
}
return(TRUE)
}
parse.AFNI.args <- function ( args, params = NULL,
other_ok=TRUE,
verb = 0) {
#for (i in 1:length(args)) {
# cat (i, args[[i]],'\n');
#}
if (is.null(args)) return(NULL)
if (!is.null(params)) {
allowed_options <- sort(names(params));
duplicate_okvec <- vector('character');
for (i in 1:1:length(params)) {
pl <- params[i][[1]];
if (pl['duplicate_ok'][[1]]) {
duplicate_okvec <- c(duplicate_okvec, names(params)[i])
}
}
} else {
allowed_options <- vector('character');
}
#Add global args like -overwrite
allowed_options <- c(allowed_options, "-overwrite")
#find locations of -*
ii <- grep ('^-.*', args);
iflg <- vector('numeric')
if (length(ii) > 0) {
for (i in 1:length(ii)) {
if (!is.num.string(args[[ii[i]]])) {
if (!length(allowed_options)) {
iflg <- append(iflg, ii[i]);
} else { #Make sure it is an acceptable name
if (length(which(args[[ii[i]]] == allowed_options))) {
iflg <- append(iflg, ii[i]);
}
}
}
}
}
if (verb) note.AFNI(paste(args[iflg]))
ops = list()
used <- vector('logical', length(args));
if (length(iflg)) {
iflg <- append(iflg,length(args)+1)
#store results
nm <- vector('character');
for (i in 1:(length(iflg)-1)) {
if (0) { # Why remove the -?, makes things inconsistent elsewhere
#newnm <- strsplit(args[[iflg[i]]],'-')[[1]][2]
} else newnm <- args[[iflg[i]]]
if (length(nm) && length(which(newnm == nm)) &&
(newnm != '-gltLabel') && (newnm != '-gltCode') && # 10/18/2012 GC: added this line for 3dMVM
(newnm != '-glfLabel') && (newnm != '-glfCode') && # 12/22/2014 GC: added this line for 3dMVM
(!length(duplicate_okvec) ||
length(which(iflg[i] == duplicate_okvec))) ){
warning(paste('option ', newnm, 'already specified.\n'),
immediate. = TRUE);
show.AFNI.args(ops)
return(NULL);
}
#nm <- append(nm, newnm)
used[iflg[i]] <- TRUE;
istrt = iflg[i]+1;
pp <- vector('character');
if (istrt <= length(args) && istrt != iflg[i+1]) {
iend <- max(c(iflg[i+1]-1, istrt))
for (ii in istrt:iend) {
pp <- append(pp, args[[ii]]);
used[ii] <- TRUE;
}
}
#create a cleaned up string
pp <- paste(pp, collapse = ' ')
if ( length(grep('^".*"$',pp)) || #Quoted string, do not split
length(grep("^'.*'$",pp)) || #Quoted string, do not split
length(grep("^1D:.*$",pp)) || #1D: string, do not split
length(grep("^R:.*$",pp)) ) {
} else {
if (verb) {
note.AFNI(sprintf("Splitting >>>%s<<<", pp))
}
pp <- strsplit(clean.args.string(pp), ' ')
}
if((newnm != '-gltLabel') && (newnm != '-glfLabel') && (newnm != '-gltCode') && (newnm != '-glfCode')) { # 12/22/2014 GC: added this line for 3dMVM
ops <- c(ops, (pp))
names(ops)[length(ops)] <- newnm
} else if(length(which(newnm == nm)))
ops[[newnm]] <- c(ops[[newnm]], pp) else {
ops <- c(ops, list(pp))
names(ops)[length(ops)] <- newnm
}
nm <- append(nm, newnm)
}
}
#cleanup
if (length(ops)) {
for (i in 1:length(ops)) {
#ops[[i]] <- clean.args.string(ops[[i]])
}
}
#numeric changes
if (length(ops))
for (i in 1:length(ops))
if(!is.list(ops[[i]])) if (is.num.string(ops[[i]]))
ops[[i]] <- as.numeric(ops[[i]])
#defaults
pp <- c(args[used == FALSE])
ops <- c (ops, list("other"=pp));
#add allowed options
ops <- c (ops, list("allowed_options"=allowed_options));
if (!other_ok) {
if (length(ops[['other']])) {
err.AFNI(paste('Illegal parameters on command line:\n',
' ', ops['other'],
'\nTry -allowed_options, or -help for details\n',
'\n'));
exit.AFNI(1);
}
}
#check
if (!check.AFNI.args(ops, params)) {
return(NULL);
} else {
return(ops);
}
}
AFNI.new.options.list <- function(history = '', parsed_args = NULL) {
lop <- list (com_history = history);
#Look for defaults
lop$overwrite <- FALSE
for (i in 1:length(parsed_args)) {
opname <- strsplit(names(parsed_args)[i],'^-')[[1]];
opname <- opname[length(opname)];
switch(opname,
overwrite = lop$overwrite <- TRUE )
}
return(lop)
}
#------------------------------------------------------------------
# Some utilities
#------------------------------------------------------------------
#history.AFNI (is almost identical to history() but without
#the interactive mode and grep
history.AFNI <- function (max.show = 25, reverse = FALSE, pattern, ...)
{
file1 <- tempfile("Rrawhist")
savehistory(file1)
rawhist <- readLines(file1)
unlink(file1)
if (!missing(pattern))
rawhist <- unique(grep(pattern, rawhist, value = TRUE,
...))
nlines <- length(rawhist)
if (nlines) {
inds <- max(1, nlines - max.show):nlines
if (reverse)
inds <- rev(inds)
}
else inds <- integer(0L)
file2 <- tempfile("hist")
writeLines(rawhist[inds], file2)
#file.show(file2, title = "R History", delete.file = TRUE)
mm <- readLines(file2)
file.remove(file2)
cat(mm, sep='\n')
invisible(mm)
}
hgrep <- function (pattern=NULL){
if (is.null(pattern)) {
history.AFNI(max.show=200)
} else {
history.AFNI(max.show=Inf, pattern=pattern)
}
}
hsgrep <- function (pattern='source'){
hgrep(pattern)
}
#Report objects using the most memory
R.bit.version <- function () {
if (.Machine$sizeof.pointer == 4) return(32)
else if (.Machine$sizeof.pointer == 8) return(64)
else return(0)
}
memory.hogs <- function (n=10, top_frac=0.9, test=FALSE, msg=NULL) {
if (test) {
toy <- function() { g <- array(0,c(128,128,128,10));
memory.hogs(test=FALSE) }
moy <- function() {k <- array(0,c(128, 128, 128, 3)); toy() }
moy()
return(NULL)
}
us <- 1024*1024
nfr <- sys.nframe()
z<- NULL
pp <- 0
while (nfr > 0) {
nfr <- nfr -1
pp <- pp+1
envir <- parent.frame(pp)
envfunc <- (as.character(sys.call(-(pp))))
if (!length(envfunc)) envfunc <- 'GLOBAL'
zs <- sapply(ls(envir), function(x)
object.size(get(x, envir = envir)))
names(zs) <- paste(names(zs)," (in ",envfunc,")")
z <- c(z, zs)
}
#Stop at the top_frac
if (top_frac > 0.0) {
m <- rev(sort(z))[1:n]
mtot <- sum(m)
ms <- m[1]
nf <- 1
while (ms /mtot < top_frac && nf <= n) {
nf <- nf+1
ms <- ms+m[nf]
}
}
m <- as.matrix(rev(sort(z))[1:nf]/us)
colnames(m) <- c("Size(Mb)")
m <- rbind(m,' total' = sum(m));
m <- rbind(m,' Grand total' = sum(z)/us);
cat( sprintf("Memory hogs check from function: %s\nR-%d bit version\n",
as.character(paste(sys.call(-1), collapse=' ')),
R.bit.version()) );
if (!is.null(msg)) cat(msg)
print(m, digits=2)
invisible(m)
}
newid.AFNI <-function(len=26) {
return(paste("XYZ_",
paste(sample(c(rep(0:9,each=5),LETTERS, letters),
len-4, replace=TRUE), collapse='', sep=''), sep='',collapse=''))
}
sys.AFNI <- function(com=NULL, fout=NULL, ferr=NULL, echo=FALSE,
CoS=errex.AFNI ) {
ss <- list(stat=1, err='', out='')
if (is.null(com)) return(ss)
rmfout <- FALSE
if (is.null(fout)) {
fout <- sprintf('/tmp/fout.%s', newid.AFNI())
rmfout <- TRUE
}
rmferr <- FALSE
if (is.null(ferr)) {
ferr <- sprintf('/tmp/ferr.%s', newid.AFNI())
rmferr <- TRUE
}
if (echo) note.AFNI(paste("Executing shell command",com))
ss$stat<-try(system(paste(com,' 1>', fout,' 2>', ferr, collapse='')))
ss$out <- readLines(fout, n=-1, ok=TRUE, warn=FALSE, encoding='unknown')
if (rmfout) unlink(fout);
ss$err <- readLines(ferr, n=-1, ok=TRUE, warn=FALSE, encoding='unknown')
if (rmfout) unlink(ferr);
if (ss$stat && !is.null(CoS)) CoS(paste("Error status executing:\n",com))
invisible(ss)
}
who.called.me <- function (quiet_inquisitor=FALSE, trim = 0) {
mm <- (as.list(sys.calls()))
#str(mm)
N_mm <- length(mm)
callstr <- NULL
if (quiet_inquisitor) skp <- 2
else skp <- 1
callstr <- ''
for (i in (N_mm-skp):1 ) {
caller <- as.character(mm[i])
if (length(caller) == 0) {
caller <- 'R_prompt'
}
if (trim==-1) { #function only
caller <- strsplit(caller[1],'(', fixed=TRUE)[[1]][1]
} else if (trim>0) {
fun <- strsplit(caller[1],'(', fixed=TRUE)[[1]][1]
par <- strsplit(caller[1],'(', fixed=TRUE)[[1]][2]
par <- strsplit(par,'')[[1]]
n_char <- min(length(par),trim)
if (n_char < length(par)) ell <- '...'
else {
ell <- ''
#and remove last parentheses
par <- sub(')$','',par)
}
if (is.na(par[1])) {
caller <- paste(fun, collapse='', sep='')
} else {
caller <- paste(fun, '(', paste(par[1:n_char],collapse='', sep=''),
ell, ')', collapse='', sep='')
}
}
spc = ' '
if (i == N_mm-skp) {
callstr <- paste(callstr, caller[1])
} else {
spc <- paste(spc, ' ', sep='')
callstr <-paste(callstr, '\n',spc, '-->',
caller,
sep= '')
}
}
#if (BATCH_MODE) caller <- as.character(sys.call(-1))
#else caller <- as.character(sys.call(-2))
return(callstr)
}
#return 1 if all strings in vector ss can be changed to numbers
is.num.string <- function(ss) {
if (is.null(ss) || !length(ss) || ss == '' ||
is.null(tryCatch(as.numeric(ss),
warning=function(ex) {}))) {
return(0);
} else {
return(1);
}
}
clean.args.string <- function(ss) {
if (is.list(ss) || length(ss) > 1) {
warning(paste('Function only works on single strings',
str(ss),'\n', sep=''),
immediate.=TRUE);
return(NULL);
}
#remove trailing whites
ss <- sub('^[[:space:]]*','',ss);
ss <- sub('[[:space:]]*$','',ss);
#remove multiple whites
ss <- gsub('[[:space:]]+',' ',ss);
#treat = nicely
ss <- gsub('[[:space:]]*=[[:space:]]*','=',ss)
return(ss)
}
deblank.string <- function(s, start=TRUE, end=TRUE, middle=FALSE) {
if (end) {
s = sub('[[:space:]]+$','',s);
}
if (start) {
s = sub('^[[:space:]]+','',s);
}
if (middle) {
s = gsub('[[:space:]]+',' ',s);
}
return(s);
}
pad.string.lines <- function(s, pre=' ', post=NULL) {
if (!is.null(pre)) {
s = sub('^',pre,s);
s = gsub('\n',sprintf('\n%s',pre),s);
}
if (!is.null(post)) {
s = sub('$',post,s);
s = gsub('\n',sprintf('%s\n',post),s);
}
return(s);
}
trim.string <- function (s, nchar=32, left=TRUE, strim='...')
{
ss <- strsplit(s,'')[[1]]
if (length(ss)>nchar) {
#try deblanking
s <- deblank.string(s)
ss <- strsplit(s,'')[[1]]
nc <- length(ss)
if (nc>nchar) {
#browser()
nstrim = length(strsplit(strim,'')[[1]])
if (left) {
ns <- nc - nchar - nstrim
if (ns > nstrim) {
ss <- ss[ns:nc]
s<-paste(strim,paste(ss,collapse=''), sep='')
}
return(s)
}else {
ns <- nchar - nstrim
ss <- ss[1:ns]
s<-paste(paste(ss,collapse=''), strim, sep='')
}
} else return(s)
} else return(s)
}
as.num.vec <- function(ss, addcount=TRUE, sepstr='.', reset=FALSE) {
if (is.list(ss) || length(ss) > 1) {
warning(paste('Function only works on single strings',
str(ss),'\n', sep=''),
immediate.=TRUE);
return(NULL);
}
ss <- clean.args.string(ss)
dd <- strsplit(ss,' ')[[1]];
nn <- vector('numeric');
ww <- vector('character');
lastname <- '.v'
valnum <- 0
for (ii in 1:length(dd)) {
vv <- strsplit(dd[ii],'=')[[1]];
if (length(vv) > 1) {
valnum <- valnum+1
ll <- vv[1]
vv <- as.numeric(vv[length(vv)]);
if (is.na(vv)) { return(NULL); }
lastname <- ll
} else {
valnum <- valnum+1
wrn <- getOption('warn'); options(warn=-1);
vv <- as.numeric(vv[1]); options(warn=wrn);
if (is.na(vv)) { return(NULL); }
if (addcount) {
sfnd <- paste(lastname, sepstr,'[[:digit:]]*$', sep='', collapse='')
if (!reset) {
ifnd <- grep(sfnd,ww);
} else {
ifnd <- grep(sfnd,ww[length(ww)]);
if (length(ifnd)) {
ifnd <- length(ww);
} else {
valnum <- 1
}
}
if (length(ifnd)) {
lastval <- strsplit(ww[ifnd[length(ifnd)]],
paste(lastname, sepstr,sep=''))[[1]];
if (lastval[length(lastval)] == '') {
valnum <- 1
} else {
valnum <- as.numeric(lastval[length(lastval)]) + 1
}
}
ll <- paste(lastname,sepstr, as.numeric(valnum), sep='')
} else {
ll <- paste(lastname, sep='')
}
}
nn <- c(nn,vv)
ww <- c(ww,ll)
}
names(nn) <- ww
return(nn)
}
as.char.vec <- function(ss) {
if (is.list(ss) || length(ss) > 1) {
warning(paste('Function only works on single strings',
str(ss),'\n', sep=''),
immediate.=TRUE);
return(NULL);
}
ss <- clean.args.string(ss)
dd <- strsplit(ss,' ')[[1]];
nn <- vector('character');
ww <- vector('character');
for (ii in 1:length(dd)) {
vv <- strsplit(dd[ii],'=')[[1]];
if (length(vv) > 1) ll <- vv[1]
else ll <- paste('v',as.character(vv[1]), sep='')
vv <- as.character(vv[length(vv)]);
if (is.na(vv)) { return(NULL); }
nn <- c(nn,vv)
ww <- c(ww,ll)
}
names(nn) <- ww
return(nn)
}
#------------------------------------------------------------------
# Functions to read 1D and other tables
#------------------------------------------------------------------
expand_1D_string <- function (t) {
vvf = vector(length = 0, mode="numeric")
#replace .. with : and split at ,
s = strsplit(sub("..",":",t, fixed=TRUE), ",")[[1]]
#Now loop and form vector of components
for (si in s) {
#cat ("working ", si, "\n")
if (length(grep(":",si))) {
vv <- eval(parse(text=si))
} else if (length(grep("@",si))) {
ssi = as.numeric(strsplit(si,"@")[[1]])
#cat ("ssi = ", ssi, "\n")
vv <- rep(ssi[2], ssi[1])
} else {
vv <- as.numeric(si)
}
#cat(si," = ",vv, "\n")
vvf <- c(vvf, vv)
#cat("vvnow = ",vvf, "\n")
}
return(vvf)
}
is.wholenumber.AFNI <- function(x, tol = .Machine$double.eps^0.5, all=TRUE) {
if (is.null(x)) return(FALSE)
if (!all) {
return(abs(x - round(x)) < tol)
} else {
return(prod(abs(x - round(x)) < tol)==1)
}
return(FALSE)
}
r.NI_new_element <- function (name, dat=NULL, atlist=NULL, tp=NULL){
nel <- list (name=name, atlist=atlist, dat=NULL)
if (is.null(tp)) {
if (is.character(dat)) {
tp <- "String"
} else if (is.numeric(dat)) {
if (is.wholenumber.AFNI(dat)) tp<-"int"
else tp <- "float"
}
}
if (!is.null(dat)) {
if (is.vector(dat)) {
nel <- r.NI_set_attribute(nel, "ni_dimen", length(dat))
nel <- r.NI_set_attribute(nel, "ni_type", tp)
nel$dat <- matrix(as.character(dat), length(dat),1)
} else if (is.matrix(dat)){
nel <- r.NI_set_attribute(nel, "ni_dimen", paste(dim(dat),collapse=','))
nel <- r.NI_set_attribute(nel, "ni_type",
paste(rep(tp,dim(dat)[2]), collapse=','))
nel$dat <- matrix(as.character(dat), dim(dat)[1], dim(dat)[2])
} else {
err.AFNI("Bad dat");
}
if (tp == "String") {
nel$dat <- apply(nel$dat, c(1,2), paste)
}
}
return(nel)
}
r.NI_get_attribute <- function (nel,name, brsel=NULL,
colwise=FALSE, is1Dstr=FALSE,
sep=' ; ', num=FALSE) {
ffs <- NULL
for (i in 1:length(nel$atlist)) {
if (!is.na(nel$atlist[[i]]$lhs) &&
nel$atlist[[i]]$lhs == name) {
ffs <- nel$atlist[[i]]$rhs
break
}
}
if (is.null(ffs)) return(NULL)
#ffs <- gsub("[;\"]","", ffs)
if (colwise) { #a per column deal, process it
#remove bad chars and split into one for each column
if (is1Dstr) {
ffs = expand_1D_string (ffs)
} else {
ffsv <- strsplit(ffs,sep)[[1]]
#some components are blanks to be killed
ffs <- ffsv[which(nchar(ffsv)!=0)]
}
if (!is.null(brsel)) {
if (max(brsel+1) > length(ffs)) {
err.AFNI(paste("Have ", length(ffs),
"attribute elements in ", paste(ffs,collapse=' '),
".\nBrick selection calls for max. of ",
max(brsel+1)," columns\n"
));
}
ffs <- ffs[brsel+1]
if (num) return(as.numeric(ffs))
else return(ffs)
} else {
#No selection, return all
if (num) return(as.numeric(ffs))
else return(ffs)
}
} else {
if (num) return(as.numeric(ffs))
else return(ffs)
}
return(NULL)
}
r.NI_set_attribute <- function (nel,name, val) {
ig <- -1
if (length(nel$atlist)) {
for (i in 1:length(nel$atlist)) {
if (!is.na(nel$atlist[[i]]$lhs) &&
nel$atlist[[i]]$lhs == name) {
ig<-i
break
}
}
}
if (ig > 0) {
nel$atlist[[ig]] <- list(lhs=deblank.string(name),
rhs=r.NI_dequotestring(val, db=TRUE))
} else {
nel$atlist <- c (nel$atlist,
list(list(lhs=deblank.string(name),
rhs=r.NI_dequotestring(val, db=TRUE))))
}
return(nel)
}
r.NI_dequotestring <- function(val, db=FALSE) {
if (db) val <- deblank.string(val)
val <- sub('^\"','', val)
val <- sub('\"$','', val)
val <- sub("^\'",'', val)
val <- sub("\'$",'', val)
return((val))
}
# A very basic parser, works only on simple elements
#of ascii headersharp niml files . Needs lots of work!
r.NI_read_element <- function (fname, HeadOnly = TRUE) {
fnp <- parse.AFNI.name(fname)
if (!file.exists(fnp$file)) {
return(NULL)
}
ff <- scan(fnp$file, what = 'character', sep = '\n', quiet=TRUE)
return(r.NI_read_str_element(ff, HeadOnly))
}
r.NI_read_str_element <- function (ff, HeadOnly = TRUE) {
nel <- list(atlist=NULL, dat=NULL)
#Remove #
ff <- gsub('^[[:space:]]*#[[:space:]]*', '',ff)
if (!length(grep('^<', ff))) {
return(NULL)
}
#strt markers
strv <- grep ('<', ff)
#stp markers
stpv <- grep ('>', ff)
#check
if (length(strv) != length(stpv)) {
err.AFNI("Have extra brackets");
return(NULL)
}
if (length(strv) < 1) {
err.AFNI("Have nothing");
return(NULL)
}
#If you have a group, jump to next element.
shead <- ff[strv[1]:stpv[1]]
if (length(grep('ni_form[[:space:]]*=[[:space:]]*\"ni_group\"', shead))
&& stpv[1]+1<length(ff)) {
#try again
return(r.NI_read_str_element(ff[stpv[1]+1:length(ff)], HeadOnly))
}
#Get only the first element, for now
for (i in 1:1) {
shead <- ff[strv[i]:stpv[i]]
#get the name
nel$name <- strsplit(shead[1],'<')[[1]][2]
nel$atlist <- vector()
#remove last > from shead
shead[length(shead)] <- sub('>$','',shead[length(shead)])
for (j in 2:length(shead)) {
attr <- strsplit(shead[j],'=')[[1]]
if (length(attr) == 1) {
nel <- r.NI_set_attribute(nel, attr[1], "")
} else if (length(attr) == 2) {
nel <- r.NI_set_attribute(nel, attr[1], attr[2])
} else if (length(attr) > 2) {
err.AFNI(paste("Parse error for ", shead[j]));
}
}
}
if (HeadOnly) return(nel)
#Get the data part, the lazy way
#All the content here is text. Not sure what to do with
#this next. Wait till context arises
tp <- r.NI_get_attribute(nel,"ni_type");
if (tp == 'String') {
for (i in (stpv[1]+1):(strv[2]-1)) {
if (is.null(nel$dat))
nel$dat <- ff[i]
else nel$dat <- c(nel$dat, ff[i])
}
} else {
for (i in (stpv[1]+1):(strv[2]-1)) {
#browser()
if (is.null(nel$dat))
nel$dat <- rbind(r.NI_dequotestring(
strsplit(deblank.string(ff[i]),split=' ')[[1]]))
else nel$dat <- rbind(nel$dat,
r.NI_dequotestring(strsplit(
deblank.string(ff[i]),split=' ')[[1]]))
}
}
return(nel)
}
r.NI_write_str_element<- function(nel, HeadOnly = TRUE) {
ff <- paste("<",nel$name,'\n', sep='')
for (i in 1:length(nel$atlist)) {
ff <- paste(ff,
' ',nel$atlist[[i]]$lhs,'="', nel$atlist[[i]]$rhs,'"\n', sep='')
}
ff <- paste(ff,'>');
if (!HeadOnly) {
if (!is.null(nel$dat)) {
ff <- paste(ff,'\n')
if (is.matrix(nel$dat)) {
for (j in 1:dim(nel$dat)[1]) {
ff <- paste(ff,paste(nel$dat[j,], collapse= ' ', sep=''),
'\n', sep='')
}
} else {
ff <- paste(ff,paste(nel$dat, collapse= '', sep=''),'\n',sep='')
}
}
}
ff <- paste(ff,'</',nel$name,'>', sep='')
return(ff)
}
is.NI.file <- function (fname, asc=TRUE, hs=TRUE) {
fnp <- parse.AFNI.name(fname)
if (!file.exists(fnp$file)) {
return(FALSE)
}
ff <- scan(fnp$file, what = 'character', nmax = 2, sep = '\n', quiet=TRUE)
if (asc && hs) {
if (!length(grep('^[[:space:]]*#[[:space:]]*<', ff))) {
return(FALSE)
} else {
return(TRUE)
}
}
return(FALSE)
}
apply.AFNI.matrix.header <- function (fname, mat,
brsel=NULL, rosel=NULL,rasel=NULL,
nheadmax = 10, checkNI=TRUE) {
attr(mat,'name') <- 'noname'
attr(mat,'FileName') <- fname
fnp <- parse.AFNI.name(fname)
if (!file.exists(fnp$file)) {
return(mat)
}
#Does this look 1D nimly?
if (checkNI && !is.NI.file(fnp$file, asc=TRUE, hs=TRUE)) {
return(mat)
}
nel <- r.NI_read_element(fnp$file, HeadOnly = TRUE)
# Common attributes
attr(mat,'name') <- nel$name
if (!is.null(labels <-
r.NI_get_attribute(nel, 'ColumnLabels', brsel, colwise=TRUE))){
if (length(labels) == ncol(mat)) colnames(mat) <- labels
}
if (nel$name == 'matrix') {
if (!is.null(colg <-
r.NI_get_attribute(nel, 'ColumnGroups', brsel,
colwise=TRUE, is1Dstr=TRUE))){
if (length(colg) == ncol(mat)) attr(mat,'ColumnGroups') <- colg
if (!is.null(labels)) {
llv <- vector(length=0,mode="numeric")
tcolg <- unique(colg)
for (i in tcolg) {
if (i>0) {
ll <- which(colg == i)
llv <- c(llv,ll[1])
}
}
tnames <-paste(strsplit(labels[llv],"#0"))
attr(mat,'TaskNames') <- tnames
}
}
if (!is.null(TR <- r.NI_get_attribute(nel, 'RowTR'))){
attr(mat, 'TR') <- as.double(TR)
}
} else if (nel$name == 'DICE') {
} else if (nel$name == '3dhistog' || nel$name == 'seg_histogram') {
if (!is.null(bw <- r.NI_get_attribute(nel, 'BinWidth'))){
attr(mat, 'BinWidth') <- as.double(bw)
}
if (!is.null(bw <- r.NI_get_attribute(nel, 'window'))){
attr(mat, 'BinWidth') <- as.double(bw)
}
if (!is.null(bw <- r.NI_get_attribute(nel, 'xlabel'))){
attr(mat, 'xlabel') <- bw
}
} else {
if (0) { #No need to whine
warn.AFNI(paste(
"Don't know what to do with attribute of element ", nel$name));
}
}
return(mat)
}
read.AFNI.xmat <- function (xmatfile, nheadmax=10) {
if (!file.exists(xmatfile)) {
err.AFNI(paste("xmatfile ",xmatfile,"not found"));
return(NULL);
}
#Now load the whole deal, comments are skipped by defaults
ff <- as.matrix(read.table(xmatfile))
#Check for attributes
ff <- apply.AFNI.matrix.header(xmatfile, ff)
#Make matrix be ts
ff <- ts(ff, names=colnames(ff), deltat=attr(ff,"TR"))
return(ff)
}
xmat.base.index <- function (xmat, AFNIindex = FALSE) {
cg <- attr(xmat,'ColumnGroups')
an <- which(cg == -1)
if (AFNIindex) an <- an -1
return(an)
}
xmat.motion.index <- function (xmat, AFNIindex = FALSE) {
cg <- attr(xmat,'ColumnGroups')
an <- which(cg == 0)
if (AFNIindex) an <- an -1
return(an)
}
xmat.roni.index <- function (xmat, AFNIindex = FALSE) {
cg <- attr(xmat,'ColumnGroups')
an <- which(cg <= 0)
if (AFNIindex) an <- an -1
return(an)
}
xmat.alltasks.index <- function (xmat, AFNIindex = FALSE) {
cg <- attr(xmat,'ColumnGroups')
an <- which(cg > 0)
if (AFNIindex) an <- an -1
return(an)
}
xmat.select.indices <- function(selstr, xmat, AFNIindex = FALSE) {
#paste(selstr, collapse=" ") is used because at times selstr is a vector
sel <- strsplit(paste(selstr, collapse=" "),
" ")[[1]]
#str(sel)
#cat(sel, 'length:' , length(sel), '\n')
ilst = vector('integer');
if (length(sel)) {
for (isel in sel) {
#cat("processing: ",isel,'\n')
if (!inherits( e <- try(as.integer(eval(parse(text=isel))),
silent=TRUE),
"try-error")) {
ilst <- append(ilst,
e+1,
after=length(ilst))
} else {
#cat("processing as string: ",isel,'\n')
if (isel == 'ALL') {
ilst <- append(ilst, seq(1,ncol(xmat)))
} else if (isel == 'ALL_TASKS') {
ilst <- append(ilst,
xmat.alltasks.index(xmat, AFNIindex=FALSE))
} else if (isel == 'RONI') {
ilst <- append(ilst, xmat.roni.index(xmat, AFNIindex=FALSE))
} else if (isel == 'MOTION') {
ilst <- append(ilst, xmat.motion.index(xmat, AFNIindex=FALSE))
} else if (isel == 'BASE') {
ilst <- append(ilst, xmat.base.index(xmat, AFNIindex=FALSE))
} else {
if (0) { #Too lose
ilst <- append(ilst, grep(isel, colnames(xmat)),
after=length(ilst))
} else {
ilst <- append(ilst,
grep(sprintf('^%s',isel), colnames(xmat)),
after=length(ilst))
}
}
}
}
} else {
ilst = 0:(ncol(xmat)-1)
}
if (length(ilst) && AFNIindex) ilst <- ilst -1
return(ilst)
}
read.AFNI.matrix.test <- function(verb=1) {
cat ( 'Running read.AFNI.matrix in test mode\n',
'See read.AFNI.matrix.test for details' )
i <- 0
mm <- paste ( 'subj age weight height\n',
'joe 13 299 123 \n',
'jane 22 600 234 \n',
'jim 2 188 23\n',
sep='', collapse='');
fouts <- sprintf('___fout%02d.1D', i);
cat (mm,file=fouts);
comm <- 'read.AFNI.matrix(fouts, verb=verb)'
note.AFNI(sprintf("Working with fouts=%s:\n %s", fouts, comm))
print(eval(parse(text=comm)))
comm <- paste("read.AFNI.matrix(fouts, verb=verb,",
"userrownames=c('jim','jane'))",
sep='', collapse='')
note.AFNI(sprintf("Working with fouts=%s:\n %s", fouts, comm))
print(eval(parse(text=comm)))
comm <- paste("read.AFNI.matrix(fouts, verb=verb, ",
"userrownames=c('jim','jane'),",
"usercolnames=c('weight','age'))",
sep='', collapse='')
print(eval(parse(text=comm)))
i<- i+1
mm <- paste ( ' age weight height\n',
' 13 299 123 \n',
' 22 600 234 \n',
' 2 188 23\n',
sep='', collapse='');
fouts <- sprintf('___fout%02d.1D', i);
cat (mm,file=fouts);
note.AFNI(sprintf("Working with %s", fouts))
comm <- 'read.AFNI.matrix(fouts, verb=verb)'
note.AFNI(sprintf("Working with fouts=%s:\n %s", fouts, comm))
print(eval(parse(text=comm)))
comm <- paste("read.AFNI.matrix(fouts, verb=verb,",
"usercolnames=c('height','weight,','age'))",
sep='', collapse='')
note.AFNI(sprintf("Working with fouts=%s:\n %s", fouts, comm))
print(eval(parse(text=comm)))
i<- i+1
mm <- paste ( ' 13 299 123 \n',
' 22 600 234 \n',
' 2 188 23\n',
sep='', collapse='');
fouts <- sprintf('___fout%02d.1D', i);
cat (mm,file=fouts);
note.AFNI(sprintf("Working with %s", fouts))
print(read.AFNI.matrix(fouts, verb=verb))
print(read.AFNI.matrix(fouts, verb=verb, usercolnames=c('col02'),
userrownames=c('row03','row01','row01')))
i<- i+1
mm <- paste ( ' heft \n',
' 299 \n',
' 600 \n',
' 188 \n',
sep='', collapse='');
fouts <- sprintf('___fout%02d.1D', i);
cat (mm,file=fouts);
note.AFNI(sprintf("Working with %s", fouts))
print(read.AFNI.matrix(fouts, verb=verb))
print(read.AFNI.matrix(fouts, verb=verb, usercolnames=c('heft')))
i<- i+1
mm <- paste ( ' 299 \n',
' 600 \n',
' 188 \n',
sep='', collapse='');
fouts <- sprintf('___fout%02d.1D', i);
note.AFNI(sprintf("Working with %s", fouts))
cat (mm,file=fouts);
print(read.AFNI.matrix(fouts, verb=verb))
i <- i+1
mm <- paste ( 'subj age \n',
'joe 13 \n',
'jane 22 \n',
'jim 2 \n',
sep='', collapse='');
fouts <- sprintf('___fout%02d.1D', i);
cat (mm,file=fouts);
comm <- 'read.AFNI.matrix(fouts, verb=verb)'
note.AFNI(sprintf("Working with fouts=%s:\n %s", fouts, comm))
print(eval(parse(text=comm)))
comm <- paste("read.AFNI.matrix(fouts, verb=verb,",
"userrownames=c('jim','jane'))",
sep='', collapse='')
note.AFNI(sprintf("Working with fouts=%s:\n %s", fouts, comm))
print(eval(parse(text=comm)))
for (j in 0:i) {
system(sprintf('\\rm -f ___fout%02d.1D', i));
}
}
AFNI.matrix.input.help.string <- function() {
s<-
'You specify 1D input in a variety of ways:
FILE.1D : Name of 1D file
1D_EXPR: A 1D expression a la AFNI "1D: 1,2,5"
R_EXPR: An R expression, "R: rep(seq(1,3),2)"
'
return(s)
}
read.AFNI.matrix <- function (fname,
usercolnames=NULL,
userrownames=NULL,
checkNI = TRUE, verb = 0) {
if (length(fname)>1) {
err.AFNI(paste("Cannot handle more than one name.\nHave ",
paste(fname,collapse=' '), sep=''))
return(NULL);
}
if (fname == '-self_test') {
read.AFNI.matrix.test()
return(NULL);
}
if (!is.null(mm <- eval.AFNI.string(fname))) {
#Might need to add some names someday...
attr(mm,'FileName') <- paste('STR<',fname,'>',sep='')
return(as.matrix(mm))
}
if (verb) cat(who.called.me())
if (is.character(fname)) {
fname <- parse.AFNI.name(fname)
fnameattr <- fname$orig_name
}else{
fnameattr <- 'NONAME'
}
if (fname$type == 'NIML') {
nmout <- sprintf('/tmp/fout.%s.1D.dset', newid.AFNI())
if (verb) warn.AFNI("Clumsy handling of NIML files via ConvertDset ...");
com <- paste(
'ConvertDset -overwrite -o_1D -prefix ', nmout, ' -i',
fname$orig_name);
ss <- sys.AFNI(com, echo = FALSE);
if (ss$stat == 1) {
err.AFNI(paste("Failed to get keys from labeltable",ltfile, ".\n",
"Command ",com,"Failed"));
return(NULL);
}
mm<-read.AFNI.matrix(nmout);
sys.AFNI(sprintf('\\rm -f %s',nmout));
return(mm);
} else {
#str(fname)
#fname$file
brk <- NULL
brk <- tryCatch({read.table(fname$file, colClasses='character')},
error=function(a){})
if (is.null(brk)) { #try as niml, just in case
if (verb) note.AFNI(paste("Attempting read as NIML."), callstr='');
fff <- r.NI_read_element(fname$file, FALSE)
if (!is.null(fff$dat)) {
brk <- as.data.frame(fff$dat, stringsAsFactors=FALSE)
checkNI = FALSE;
} else {
err.AFNI(paste("Failed to read matrix from ", fname$file,".\n"));
return(NULL);
}
}
}
if ( tolower(brk$V1[1]) == 'name' ||
tolower(brk$V1[1]) == 'subj' ||
tolower(brk$V1[1]) == '#file') {
subjCol <- brk$V1[2:dim(brk)[1]];
covNames <- paste(brk[1,2:dim(brk)[2]]);
if (dim(brk)[2]-1 == 1) {
covMatrix <- as.matrix(as.numeric(brk[2:dim(brk)[1],2]))
} else {
for (ii in 1:(dim(brk)[2]-1)) { #Add one column at a time
if (ii==1) {
covMatrix <- cbind(
as.numeric(brk[2:dim(brk)[1],2:dim(brk)[2]][[ii]]));
} else {
covMatrix <- cbind(covMatrix,
as.numeric(brk[2:dim(brk)[1],2:dim(brk)[2]][[ii]]));
}
}
}
} else {
flg <- tryCatch({as.numeric(brk$V1[1])}, warning=function(aa) {});
if (is.null(flg)) { #Just labels
covNames <- paste(brk[1,1:dim(brk)[2]]);
subjCol <- paste('row',sprintf('%02d',c(1:(dim(brk)[1]-1))), sep='')
istrt<- 2
}else { #completely naked
ccc <- trim.string(fname$prefix, nchar=10, left=FALSE, strim='..');
covNames <- paste(ccc,sprintf(' c%02d',c(1:dim(brk)[2])),sep='');
subjCol <- paste(ccc,sprintf(' r%02d',c(1:dim(brk)[1])), sep='')
istrt<- 1
}
if (dim(brk)[2] == 1) {
covMatrix <- as.matrix(as.numeric(brk[istrt:dim(brk)[1],1]))
} else {
for (ii in 1:(dim(brk)[2])) { #Add one column at a time
ccc <- tryCatch(
{as.numeric(brk[istrt:dim(brk)[1],1:dim(brk)[2]][[ii]])},
warning=function(aa) {} )
if (is.null(ccc)) {
warn.AFNI(paste("Failed to process column ",
ii-1, " in ", fnameattr,
". Using NA instead"))
ccc <- NA*vector(length=dim(brk)[1]-istrt+1)
}
if (ii==1) {
covMatrix <- cbind(ccc);
} else {
covMatrix <- cbind(covMatrix,ccc);
}
}
}
}
if (verb>2) {
note.AFNI("Browser here, not active");
#browser()
}
rownames(covMatrix) <- subjCol;
colnames(covMatrix) <- covNames;
#And now apply selectors
if (!is.null(fname$brsel)) {
sbsel <- eval.AFNI.1D.string(fname$brsel, nmax=dim(covMatrix)[2]-1)
if (min(sbsel) < 0 || max(sbsel)>=dim(covMatrix)[2]) {
err.AFNI(
sprintf('column selection outside possible range of <0,%d> in %s',
dim(covMatrix)[2]-1, fname$file));
return(NULL);
}
covMatrix <- covMatrix[,sbsel+1, drop=FALSE]
} else sbsel <- NULL
if (!is.null(fname$rosel)) {
rosel <- eval.AFNI.1D.string(fname$rosel, nmax=dim(covMatrix)[1]-1)
if (min(rosel) < 0 || max(rosel)>=dim(covMatrix)[1]) {
err.AFNI(
sprintf('row selection outside possible range of <0,%d> in %s',
dim(covMatrix)[1]-1, fname$file));
return(NULL);
}
covMatrix <- covMatrix[rosel+1,, drop=FALSE]
} else rosel <- NULL
if (!is.null(fname$rasel)) {
err.AFNI('Not ready to deal with range selection');
return(NULL);
} else rasel <- NULL
#Now, reorder per user*names
if (!is.null(userrownames)) {
dd <- userrownames[!(userrownames %in% subjCol)]
if (length(dd)) {
warning (paste('Row(s) "', paste(dd,collapse=' '),
'" do(es) not have an entry.\n'),
immediate.=TRUE);
return(NULL);
}
for (ii in 1:1:length(userrownames)) {
if (ii==1) {
mm <- rbind(covMatrix[userrownames[ii],]);
} else {
mm <- rbind(mm, covMatrix[userrownames[ii],]);
}
}
rownames(mm) <- userrownames
if(length(covNames)) colnames(mm) <- covNames
covMatrix <- mm
}
if (!is.null(usercolnames)) {
dd <- usercolnames[!(usercolnames %in% covNames)]
if (length(dd)) {
warning (paste('Column(s) "', paste(dd,collapse=' '),
'" do(es) not have an entry.\n'),
immediate.=TRUE);
return(NULL);
}
for (ii in 1:1:length(usercolnames)) {
if (ii==1) {
mm <- cbind(covMatrix[,usercolnames[ii]]);
} else {
mm <- cbind(mm, covMatrix[,usercolnames[ii]]);
}
}
colnames(mm) <- usercolnames
covMatrix <- mm
}
covMatrix <- apply.AFNI.matrix.header(fnameattr, covMatrix,
brsel=sbsel, rosel=rosel,
rasel=rasel,checkNI = checkNI)
return(covMatrix)
}
read.AFNI.labeltable <- function (ltfile=NULL,
def.grp.label = c('CSF','GM','WM','InSk','OutSk','Out')) {
if (!file.exists(ltfile)) {
err.AFNI(paste(ltfile, "does not exist"));
return(NULL)
}
if (is.null(lt <- r.NI_read_element(ltfile, HeadOnly=FALSE))) {
err.AFNI("Failed to read label table");
return(NULL)
}
if (lt$name != "VALUE_LABEL_DTABLE") {
err.AFNI("lt file not VALUE_LABEL_DTABLE")
return(NULL)
}
lt$labels <- lt$dat[,2]
lt$keys <- as.numeric(lt$dat[,1])
if (is.null(gl <- r.NI_get_attribute(lt, "grp.label", colwise=TRUE))) {
gl <- def.grp.label
}
if (!is.null(gl)) {
lt$grp.label <- gl
}
grps <- r.NI_get_attribute(lt, "grp.label", colwise=TRUE, num=TRUE)
if (is.null(grps)) {
for (lab in lt$labels) {
nn <- 0
ig <- 0
for (i in 1:length(gl)) {
if (length(grep(gl[i],lab))) {
if (ig == 0) {
nn <- 1
ig <- i
} else {
if ( length(gl[ig]) != length(lab) ) {
if (length(gl[i]) == length(lab) ) {
# A better match
nn <- 1
ig <- i
} else {
#Another partial
nn <- nn + 1
ig <- i
}
}
}
}
}
if (nn != 1) {
warn.AFNI(paste(
"Failed to find 1 group for ",lab, ". Found ",
nn ,"candidates. Setting grp to",
length(gl)+1))
ig <- length(gl)+1
}
grps <- c(grps, ig)
}
}
if (!is.null(grps)) {
lt$groups <- grps
}
return(lt)
}
write.AFNI.matrix <- function (m, fname='test.1D') {
if (class(m) == "AFNI_c_dataset") {
m <- m$brk
}
if (is.vector(m)) {
write(m, fname, ncolumns=1)
} else if (is.matrix(m)) { #Need to transpose matrix
write(t(m), fname, ncolumns=dim(m)[2])
} else {
err.AFNI("Don't know how to write this thing");
}
return(fname)
}
dset.nonzero.indices <- function (dset) {
nn <- dset.ndim(dset)
if (nn == 3 || nn==1) {
return(which(rowSums(abs(dset$brk), dims=nn)!=0,arr.ind=TRUE))
} else {
err.AFNI("Bad dim");
return(NULL);
}
}
dset.zero.indices <- function (dset) {
nn <- dset.ndim(dset)
if (nn == 3 || nn==1) {
return(which(rowSums(abs(dset$brk), dims=nn)==0,arr.ind=TRUE))
} else {
err.AFNI("Bad dim");
return(NULL);
}
}
dset.nvox <- function (dset) {
nvox = -1
if (length(dd <- dim(dset$brk)) <= 2) {
nvox = dd[1]
} else if (length(dd) > 2) {
nvox = prod(dd[1:3])
}
return(nvox)
}
dset.nval <- function (dset) {
nval = -1
if (length(dd <- dim(dset$brk)) == 1) {
nval = 1
} else if (length(dd) == 2) {
nval = dd[2]
} else if (length(dd) == 3) {
nval = 1
} else if (length(dd) == 4) {
nval = dd[4]
}
return(nval)
}
dset.nx <- function(dset) {
return(dset$dim[1])
}
dset.ny <- function(dset) {
return(dset$dim[2])
}
dset.nz <- function(dset) {
return(dset$dim[3])
}
dset.nslice <- function(dset) {
return(dset$dim[3])
}
dset.ndim <- function (dset) {
ndims <- -1
if (length(dd <- dim(dset$brk)) == 1 ||
length(dd) == 2) {
ndims <- 1
} else if (length(dd) == 3 || length(dd) == 4) {
ndims <- 3
}
return(ndims)
}
get.c.AFNI.attribute<- function (atlist, attribute, asis = FALSE, verb=0) {
if (!is.null(atlist$NI_head)) atlist <- atlist$NI_head
nel <- atlist[[attribute]]
if (!is.null(nel)) {
if (verb > 1) note.AFNI(c("Found", attribute, '\n'))
if (asis) { return(nel$dat) }
else {
tp <- r.NI_get_attribute(nel,"ni_type");
if (tp == 'String') {
dd <- sub("^ *\"",'', nel$dat)
dd <- sub(" *\"$",'', dd)
return(strsplit(dd,'~')[[1]])
} else if (tp == 'float') {
return(as.numeric(nel$dat));
} else if (tp == 'int') {
return(as.integer(nel$dat));
} else {
cat ("What of ", tp)
#browser()
nel
return(NULL);
}
}
}
if (verb) {
err.AFNI(paste("Attribute:",attribute,"Not found"))
show.dset.attr(atlist)
}
return(NULL);
}
set.c.AFNI.attribute<- function (atlist, attribute, val=NULL, verb=0, tp=NULL,
strsep="~") {
if (!is.null(atlist$NI_head)) atlist <- atlist$NI_head
nel <- atlist[[attribute]]
if (!is.null(nel)) {
if (verb) cat ("Found", attribute, '\n')
tp <- r.NI_get_attribute(nel,"ni_type");
if (tp == 'String') {
nel$dat <- paste("\"",paste(val, collapse=strsep), "\"", sep='')
nel <- r.NI_set_attribute(nel, "ni_dimen", 1);
} else if (tp == 'float' || tp == 'int') {
ff <- NULL
if (is.matrix(val)) {
nel <- r.NI_set_attribute(nel,
"ni_dimen",paste(dim(val),collapse=','))
nel$dat <- matrix(as.character(val), dim(val)[1], dim(val)[2])
} else if (is.vector(val)) {
nel <- r.NI_set_attribute(nel,"ni_dimen",paste(length(val)))
nel$dat <- matrix(as.character(val), length(val),1)
} else {
nel$dat <- val
}
} else {
cat ("What of ", tp)
#browser()
nel
}
} else {
nel <- r.NI_new_element("AFNI_atr", dat=val, atlist= NULL, tp=tp)
nel <- r.NI_set_attribute(nel, "atr_name", attribute)
atlist[[attribute]] <- nel
#Now that name and type are set, call again because
#you want the multi-col strings to be turned into one
#part with the '~' delimiters....
tp <- r.NI_get_attribute(nel,"ni_type");
if (tp == 'String') {
atlist <- set.c.AFNI.attribute(atlist, attribute, val=val,
verb=verb,strsep=strsep)
return(atlist);
}
}
atlist[[attribute]] <- nel
if (verb && is.null(nel)) {
err.AFNI(paste("Attribute:",attribute,"Not found."))
show.dset.attr(atlist)
}
return(atlist);
}
#See README.attributes's SCENE_DATA and niml_stat.c's distname[]
statsym.distold2new <- function(dist) {
if (dist == 'fico') return("Correl")
if (dist == 'fitt') return("Ttest")
if (dist == 'fift') return("Ftest")
if (dist == 'fizt') return("Zscore")
if (dist == 'fict') return("Chisq")
if (dist == 'fibt') return("Beta")
if (dist == 'fibn') return("Binom")
if (dist == 'figt') return("Gamma")
if (dist == 'fipt') return("Poisson")
return(dist)
}
statsym.list2code <- function(statsym) {
if (is.null(statsym) || length(statsym)==0) return("")
imx <- statsym[[1]]$sb
if(length(statsym)>1) for (i in 2:length(statsym)) {
if (imx < statsym[[i]]$sb) imx <- statsym[[i]]$sb
}
code <- NULL
for (i in 1:imx+1) {
code <- c(code,"none");
}
for (i in 1:length(statsym)) {
code[statsym[[i]]$sb+1] <-
paste(statsym.distold2new(statsym[[i]]$typ),"(",
paste(statsym[[i]]$par,collapse=','), ")",
sep = '');
}
return(code)
}
header.version <- function(hh, defmeth='UNKNOWN') {
if (defmeth == 'AUTO') {
if (have_R_io()) defmeth <- 'clib' else defmeth <- 'Rlib'
}
if (!is.null(hh) && length(hh)>0) {
if (is.null(names(hh[[1]]))) return('Rlib');
if (length(which("atlist" == names(hh[[1]]))) > 0) return('clib');
#Cannot tell, return default
return(defmeth)
}
return(defmeth)
}
#This function is to be used in different ways depending on whether you
#are setting, or retrieving attributes
#Example:
#rs <- read.c.AFNI('littleone.niml.dset')
#to retrieve attributes:
# dset.attr(rs,"IJK_TO_DICOM")
# or
# dset.attr(rs,"BRICK_LABS")
#To set attributes:
# rs <- dset.attr(rs,"BRICK_LABS",val=c('one', 'two', 'three'))
#To set an attribute to a default if none exists
# rs <- dset.attr(rs,"BRICK_LABS",default=c('set if none existed'))
#In query mode only, one can also pass the old '$header' list
#generated by the old read.AFNI
dset.attr <- function (dset, name=NULL, colwise=FALSE, num=FALSE,
val=NULL, default=NULL, tp=NULL) {
attr <- NULL
if (class(dset) != 'AFNI_R_dataset' && class(dset) != 'AFNI_c_dataset') {
if (header.version(dset) == 'Rlib') {
#OK to change nature of dset, only queries allowed for old style
dset <- list(header=dset);
}
}
if (!is.null(dset$header)) { #olde way
if (!is.null(val)) {
err.AFNI("Cannot set values with old header format");
return(attr)
}
iattr <- which(name == names(dset$header))
if (!length(iattr)) {
return(attr)
} else {
attr <- dset$header[[iattr]]
}
if (colwise) {
attr <- sub("^'",'', attr)
attr <- strsplit(attr,'~')[[1]]
}
if (num) attr <- as.numeric(attr)
} else {
if (!is.null(dset$NI_head)) hatr <- dset$NI_head
else hatr <- dset
if (is.null(name)) {
return(names(hatr))
}
if (is.null(val) && is.null(default)) { #Retrieval
if (name == 'dim') {
dd <- get.c.AFNI.attribute(hatr, "DATASET_DIMENSIONS")
ee <- get.c.AFNI.attribute(hatr, "DATASET_RANK")
return(c(dd[1], dd[2], dd[3], ee[2]))
} else if (name == 'orient') {
return(orcode.AFNI(get.c.AFNI.attribute(hatr, "ORIENT_SPECIFIC")))
} else if (name == 'origin') {
return(get.c.AFNI.attribute(hatr, "ORIGIN"))
} else if (name == 'delta') {
return(get.c.AFNI.attribute(hatr, "DELTA"))
} else if (name == 'hist') {
return(get.c.AFNI.attribute(hatr, "HISTORY_NOTE"))
} else if (name == "TR") {
dd <- get.c.AFNI.attribute(hatr, "TAXIS_NUMS")
ee <- get.c.AFNI.attribute(hatr, "TAXIS_FLOATS")
if (dd[3] == 77001) ee[2] = ee[2]/1000; #From msec to sec
return(ee[2])
} else {
return(get.c.AFNI.attribute(hatr, name))
}
} else if (!is.null(val)) { #Set the attribute
if (name == 'note') {
i <- get.c.AFNI.attribute(hatr, "NOTES_COUNT")
if (is.null(i)) i<-0
i<-i+1
hatr <- set.c.AFNI.attribute(hatr, "NOTES_COUNT", val = i)
name <- sprintf('NOTE_NUMBER_%03d', i)
hatr <- set.c.AFNI.attribute(hatr, name, val)
} else if (name == 'dim') {
if (is.null(val)) {
if (is.null(dset$brk)) {
err.AFNI("Cannot set dims from nothing");
return(NULL);
}
val <- dset.dimBRKarray(dset)
}
if (length(val) != 4) {
err.AFNI("Bad val for dim");
return(NULL);
}
hatr <- set.c.AFNI.attribute(hatr, "DATASET_DIMENSIONS",
val=c(val[1:3], 0, 0))
hatr <- set.c.AFNI.attribute(hatr, "DATASET_RANK",
val=c(3, val[4], rep(0,6)))
} else if (name == "statsym") {
hatr <- set.c.AFNI.attribute(hatr, "BRICK_STATSYM",
statsym.list2code(statsym=val), strsep=';')
} else if (name == "TR") {
dd <- get.c.AFNI.attribute(hatr, "TAXIS_NUMS");
if (is.null(dd)) {
rr <- get.c.AFNI.attribute(hatr, "DATASET_RANK")
if (is.null(rr)) {
err.AFNI("Cannot set TAXIS_NUMS with DATASET_RANK missing");
return(NULL);
}
dd <- c(rr[2], 0, 77002, -999, -999, -999, -999, -999);
}
ee <- get.c.AFNI.attribute(hatr, "TAXIS_FLOATS");
if (is.null(ee)) {
ee <- c(0, 0, 0, 0, 0, -999999.000, -999999.000, -999999.0000);
}
dd[2] <- 77002; #seconds
ee[2] <- val;
hatr <- set.c.AFNI.attribute(hatr, "TAXIS_NUMS", val=dd);
hatr <- set.c.AFNI.attribute(hatr, "TAXIS_FLOATS", val = ee);
} else{
hatr <- set.c.AFNI.attribute(hatr, name, val, tp=tp)
}
if (!is.null(dset$NI_head)) {
dset$NI_head <- hatr;
return(dset);
} else return(hatr);
} else if (!is.null(default)) {
if (is.null(get.c.AFNI.attribute(hatr,name))) {
hatr <- set.c.AFNI.attribute(hatr, name, val=default, tp=tp)
}
if (!is.null(dset$NI_head)) {
dset$NI_head <- hatr;
return(dset);
} else return(hatr);
}
}
return(attr)
}
show.dset.attr <- function (atlist, val=FALSE, sel=NULL ) {
if (!is.null(atlist$NI_head)) atlist <- atlist$NI_head
j = 1;
while (j<=length(atlist)) {
nm <- r.NI_get_attribute(atlist[[j]],"atr_name")
if (is.null(sel) || length(grep(sel,nm, ignore.case=TRUE))) {
cat(nm)
if (val) {
cat('=', atlist[[j]]$dat,'\n')
} else {
cat('\n');
}
}
j <- j + 1;
}
invisible(NULL);
}
dset.labels <- function (dset) {
return(dset.attr(dset, 'BRICK_LABS', colwise=TRUE))
}
dimBRKarray <- function(brk=NULL) {
d <- c(1,1,1,1)
if (is.array(brk)) dd <- dim(brk)
else if (is.vector(brk)) dd <- length(brk)
else {
note.AFNI("Don't know what to do with this brk")
return(NULL)
}
d[1:length(dd)]<-dd
return(d)
}
subBRKarray <- function(brk=NULL, sel=NULL) {
d <- dimBRKarray(brk)
if (is.null(d)) {
note.AFNI("Failed to get dims")
return(NULL)
}
if (is.null(sel) || d[4] == 1) {
v <- brk
} else {
if (max(sel) > d[4] || min(sel) < 1) {
note.AFNI(printf("sel outside of range 1:%d\n", d[4]))
return(NULL)
}
v <- brk[,,,sel]
d[4] <- length(sel)
dim(v) <- d
}
return(v)
}
dset.dimBRKarray <- function(dset) {
return(dimBRKarray(dset$brk))
}
index3Dto1D <- function(ii, dd) {
dd12 <- dd[1]*dd[2]
return( (ii[,1]-1)+
(ii[,2]-1)*dd[1]+
(ii[,3]-1)*dd12 + 1)
}
index1Dto3D <- function(ii, dd) {
dd12 <- dd[1]*dd[2]
ii <- ii-1
k <- ii %/% dd12
j <- ii %% dd12
i <- j %% dd[1]
j <- j %/% dd[1]
return(cbind(i+1,j+1,k+1))
}
dset.index3Dto1D <- function(ii, dset) {
dd <- dset$dim
return( index3Dto1D(ii,dd ))
}
dset.index1Dto3D <- function(ii, dset) {
dd <- dset$dim
return( index1Dto3D(ii,dd ))
}
dset.3DBRKarrayto1D <- function(dset) {
if (is.null(dset$brk)) {
err.AFNI("NULL or non existent dset$brk")
return(NULL);
}
if (dset.ndim(dset) == 1) return(dset)
dd <- dset.dimBRKarray(dset)
dim(dset$brk) <- c(prod(dd[1:3]),dd[4])
return(dset)
}
dset.1DBRKarrayto3D <- function(dset, dd=NULL) {
if (is.null(dset$brk)) {
err.AFNI("NULL or non existent dset$brk")
return(NULL);
}
if (dset.ndim(dset) == 3) return(dset)
ddi <- dset.dimBRKarray(dset)
if (is.null(dd)) { #get it from voxel header
dd <- dset$dim
}
if (prod(dd[1:3]) != prod(ddi[1])) {
err.AFNI("Dimension mismatch");
return(NULL)
}
dim(dset$brk) <- c(dd[1:3],ddi[2])
return(dset)
}
array2dset <- function (brk=NULL, format='BRIK', meth='AUTO') {
if (meth == 'AUTO') {
if (have_R_io()) meth <- 'clib' else meth <- 'Rlib'
}
if (is.vector(brk)) {
brk <- as.array(brk, dim=c(length(brk),1,1,1))
}
dd <- dimBRKarray(brk)
if (meth=='Rlib') {
z <- list(brk=brk,format=format, dim=dd,
delta=NULL, origin=NULL, orient=NULL)
class(z) <- "AFNI_R_dataset"
} else if (meth == 'clib') {
z <- list(brk=brk, NI_head<- list())
z$NI_head <- dset.attr(z$NI_head,
"IDCODE_STRING", val="SET_AT_WRITE_RANDOM")
z$NI_head <- dset.attr(z$NI_head, "ORIGIN", val=c(0,0,0), tp="float")
z$NI_head <- dset.attr(z$NI_head, "DELTA", val=c(1,1,1), tp="float")
z$NI_head <- dset.attr(z$NI_head, "ORIENT_SPECIFIC",
val = orcode.AFNI('RAI'))
if (!is.null(names(brk)))
z$NI_head <- dset.attr(z$NI_head, "BRICK_LABS",
val = names(brk));
class(z) <- "AFNI_c_dataset"
} else {
err.AFNI('bad meth');
return(NULL)
}
return(z)
}
typecode.AFNI <- function(typestr, def='MRI_short') {
if (is.null(typestr)) return(typecode.AFNI(def))
if (is.character(typestr)) {
if (typestr == 'MRI_byte' || typestr == 'byte') return(0)
if (typestr == 'MRI_short' || typestr == 'short') return(1)
if (typestr == 'MRI_int' || typestr == 'int') return(2)
if (typestr == 'MRI_float' || typestr == 'float') return(3)
if (typestr == 'MRI_double' || typestr == 'double') return(4)
if (typestr == 'MRI_complex' || typestr == 'complex') return(5)
if (typestr == 'MRI_rgb' || typestr == 'rgb') return(6)
err.AFNI(paste('Bad typecode ', typestr));
return(NULL);
} else {
if (typestr < 0 || typestr >6) {
err.AFNI(paste('Bad typecode ', typestr));
return(NULL);
}
return(typestr);
}
err.AFNI(paste('Should not be here ', typestr));
return(NULL);
}
orcode.AFNI <- function(orstr) {
if (is.character(orstr)) {
orcode <- c (-1,-1,-1)
orstr <- strsplit(orstr,'')[[1]]
for (i in 1:3) {
switch (tolower(orstr[i]),
r = orcode[i] <- 0,
l = orcode[i] <- 1,
p = orcode[i] <- 2,
a = orcode[i] <- 3,
i = orcode[i] <- 4,
s = orcode[i] <- 5)
if (orcode[i] < 0) {
err.AFNI(paste('Bad orientation code ', orstr));
return(NULL)
}
}
} else {
for (i in 1:3) {
orcode <- orstr
if (length(orcode) != 3 || orcode[i] < 0 || orcode[i] > 5) {
err.AFNI(paste('Bad orientation code ', orstr));
return(NULL)
}
}
}
return(orcode);
}
#------------------------------------------------------------------
# Extracted (and modified) from fmri library by Karsten Tabelow,
# tabelow@wias-berlin.de and Joerg Polzehl (polzehl@wias-berlin.de)
#
# Updates by ZSS & GC
#------------------------------------------------------------------
read.AFNI <- function(filename, verb = 0, ApplyScale = 1, PercMask=0.0,
meth = 'AUTO', forcedset = FALSE) {
if (meth == 'AUTO') {
if (have_R_io()) meth <- 'clib' else meth <- 'Rlib'
}
an <- parse.AFNI.name(filename);
if (verb > 1) {
show.AFNI.name(an);
}
if ( !forcedset && ((an$type == "1D" && an$ext != ".1D.dset")
|| an$type == "Rs" || an$type == "1Ds")) {
brk <- read.AFNI.matrix(an$orig_name)
z <- array2dset(brk, format=an$type)
return(z)
}
if (meth == 'clib' || an$type == 'NIML' ||
(an$type == "1D" && an$ext == ".1D.dset")) {
return(read.c.AFNI(filename, verb = verb, ApplyScale = 1, PercMask=0.0))
}
#If you have any selectors, use 3dbucket to get what you want, then read
#temp dset. This is an ugly fix for now, but will change it later if
#I/O is issue
if (verb) { cat ('Need tmp?\n'); }
rmtmp <- 0;
if (!is.null(an$brsel) || !is.null(an$rosel) || !is.null(an$rasel)) {
rmtmp <- 1;
#Changed below from >& /dev/null to > /dev/null 2>&1
#because system uses sh not tcsh
#Tx to G. Pagnoni & Co.
com <- paste ('3dcalc -overwrite -prefix ___R.read.AFNI.' ,
basename(an$pprefix),
' -a "', filename,'" -expr "a" > /dev/null 2>&1',
sep = '');
if (try(system(com)) != 0) {
warning(paste("Failed to execute:\n ", com),
immediate. = TRUE);
return(NULL);
}
an$pprefix <- paste('___R.read.AFNI.',basename(an$pprefix), sep = '');
if (!(exists.AFNI.name(head.AFNI.name(an)))) {
warning(paste("Failed to create: ",
head.AFNI.name(an), brik.AFNI.name(an), '\n'),
immediate. = TRUE);
return(NULL);
}
}
if (verb) { cat ('Checking existence\n'); }
if (!(exists.AFNI.name(head.AFNI.name(an)))) {
err.AFNI(paste("Failed to read: ", head.AFNI.name(an),
brik.AFNI.name(an)));
return(NULL);
}
#Cannot read compressed stuff (see size usage below)
if (verb) { cat ('Uncompressing\n'); }
uncompress.AFNI(head.AFNI.name(an));
conhead <- file(head.AFNI.name(an),"r")
header <- readLines(conhead)
close(conhead)
types <- NULL
args <- NULL
counts <- NULL
values <- NULL
for (i in 1:length(header)) {
if (regexpr("^type *= *", header[i]) != -1) {
tmptype <- strsplit(header[i]," *= *")[[1]][2]
types <- c(types,tmptype)
args <- c(args,strsplit(header[i+1]," *= *")[[1]][2])
tmpcounts <- as.numeric(strsplit(header[i+2]," *= *")[[1]][2])
counts <- c(counts,tmpcounts)
i <- i+3
tmpvalue <- ""
while ((regexpr("^$", header[i]) == -1) && (i <= length(header))) {
tmpvalue <- paste(tmpvalue,header[i])
i <- i+1
}
tmpvalue <- sub("^ +","",tmpvalue)
if ((tmptype == "integer-attribute") || (tmptype == "float-attribute")) {
tmpvalue <- as.numeric(strsplit(tmpvalue," +")[[1]])
}
values <- c(values,list(value=tmpvalue))
}
}
names(values) <- args
dx <- values$DATASET_DIMENSIONS[1]
dy <- values$DATASET_DIMENSIONS[2]
dz <- values$DATASET_DIMENSIONS[3]
dt <- values$DATASET_RANK[2]
scale <- values$BRICK_FLOAT_FACS
size <- file.info(brik.AFNI.name(an))$size/(dx*dy*dz*dt)
if (is.na(size)) {
err.AFNI("Failed to determine file size");
return(NULL);
}
if (regexpr("MSB",values$BYTEORDER_STRING[1]) != -1) {
endian <- "big"
} else {
endian <- "little"
}
if (min(abs(values$DELTA)) != 0) {
weights <-
abs(values$DELTA/min(abs(values$DELTA)))
} else {
weights <- NULL
}
# browser()
if (verb) { cat ('Reading Bin\n'); }
if (as.integer(size) == size) {
conbrik <- file(brik.AFNI.name(an),"rb")
# modified below by GC 12/2/2008
if (all(values$BRICK_TYPES==0) | all(values$BRICK_TYPES==1)) {
mybrk<- readBin( conbrik, "int", n=dx*dy*dz*dt, size=size,
signed=TRUE, endian=endian)
} else if (all(values$BRICK_TYPES==3)) {
mybrk<- readBin(conbrik, "numeric", n=dx*dy*dz*dt, size=size,
signed=TRUE, endian=endian) # float
} else {
err.AFNI("Cannot read datasets of multiple data types");
close(conbrik)
return(NULL);
}
close(conbrik)
dim(mybrk) <- c(dx,dy,dz,dt)
if (ApplyScale) {
if (verb) { cat ('Scaling\n'); }
#After this operation, size of mytt doubles if initially read as int
for (k in 1:dt) if (scale[k] != 0) mybrk[,,,k] <- scale[k] * mybrk[,,,k]
}
mask=NULL;
if (PercMask > 0.0) { #ZSS: Dunno what that is for.
# 0.75 was default for PercMask
mask <- array(TRUE,c(dx,dy,dz))
mask[mybrk[,,,1] < quantile(mybrk[,,,1],PercMask)] <- FALSE
}
z <- list(brk=mybrk,format=an$type,delta=values$DELTA,
origin=values$ORIGIN,
orient=values$ORIENT_SPECIFIC,
dim=c(dx,dy,dz,dt),weights=weights, header=values,mask=mask)
} else {
warning("Error reading file: Could not detect size per voxel\n")
z <- list(brk=NULL,format=an$type,delta=NULL,
origin=NULL,orient=NULL,dim=NULL,weights=NULL,
header=values,mask=NULL)
}
class(z) <- "AFNI_R_dataset"
attr(z,"file") <- paste(filename,"BRIK",sep="")
if (rmtmp == 1) {
if (verb) {
cat ('ZSS: Will remove tmp files\n');
}
#Changed below from >& to 2&>1 because system uses sh not tcsh
#Tx to G. Pagnoni & Co.
system('\\rm -f ___R.read.AFNI.* > /dev/null 2>&1');
}else{
if (verb) {
cat ('ZSS: No temps to remove\n');
}
}
invisible(z);
}
read.c.AFNI <- function(filename, verb = 0, ApplyScale = 1, PercMask=0.0) {
if (!is.loaded('R_THD_load_dset')) {
err.AFNI("Missing R_io.so");
return(NULL);
}
if (ApplyScale != 1 || PercMask != 0.0) {
err.AFNI(paste("This function is not ready to abide by ApplyScale != 1 or",
"PercMask != 0.0"));
return(NULL);
}
uopts = list( debug=verb );
rs2 <- .Call("R_THD_load_dset", name = as.character(filename), opts = uopts )
if (is.null(rs2)) return(NULL);
hatr <- parse.c.AFNI.head(rs2$head);
rs2$head <- NULL
ddb <- dset.attr(hatr, "dim")
dim(rs2$brk) <- ddb
rs2[['format']] <- 1
rs2[['delta']] <- dset.attr(hatr, "delta")
rs2[['origin']] <- dset.attr(hatr, "origin")
rs2[['orient']] <- dset.attr(hatr, "orient")
rs2[['dim']] <- dset.attr(hatr, "dim")
rs2[['NI_head']] <- hatr
rs2[['header']] <- NULL
rs2[['mask']] <- NULL
rs2[['weights']] <- NULL
class(rs2) <- "AFNI_c_dataset"
invisible(rs2)
}
parse.c.AFNI.head <- function (head) {
j <- 1
nel <- NULL;
nms <- NULL;
while (j <= length(head)) {
nelc <- r.NI_read_str_element(strsplit(head[j], "\n")[[1]], FALSE)
if (!is.null(nelc)) {
nn <- r.NI_get_attribute(nelc, "atr_name");
if (!is.null(nn)) {
nms <- c(nms, nn)
nel <- c(nel, list(nelc));
} else {
warn.AFNI("Have element that has no atr_name");
}
}
j <- j + 1;
}
if (!is.null(nel)) {
names(nel) <- nms
}
return(nel)
}
unparse.c.AFNI.head <- function (nel) {
j <- 1
head <- NULL;
nms <- NULL;
while (j <= length(nel)) {
#r.NI_write_str_element(nel[[j]],FALSE)
nm <- r.NI_get_attribute(nel[[j]],"atr_name")
if (!is.null(nm)) nms <- c(nms,nm)
else nms <- c(nms, 'noname');
head <- c(head, r.NI_write_str_element(nel[[j]],FALSE));
j <- j + 1;
}
if (length(nms)) names(head)<-nms
return(head)
}
#A funtion to create an AFNI header string
AFNIheaderpart <- function(type, name, value) {
a <- "\n"
a <- paste(a, "type = ", type, "\n", sep="")
a <- paste(a, "name = ", name, "\n", sep="")
if (regexpr("string",type) == 1) {
value <- paste("'", value, "~", sep="")
a <- paste(a, "count = ", nchar(value) - 1, "\n", sep ="")
a <- paste(a, value, "\n", sep="")
} else {
a <- paste(a, "count = ", length(value), "\n", sep ="")
j <- 0
while (j<length(value)) {
left <- length(value) - j
if (left>4) left <- 5
a <- paste(a, paste(value[(j+1):(j+left)],collapse=" "), "\n", sep=" ")
j <- j+5
}
}
return(a);
}
#Calculate the min, and max of BRIK data y
minmax <- function(y) {
r <- NULL;
for (k in 1:dim(y)[4]) {
r <- c(r,min(y[,,,k]),max(y[,,,k]))
};
return(r);
}
newid.AFNI.old <- function(ext=0) {
if (ext) { #in house
return(
paste('GCR_',paste(
sample(c(rep(0:9,each=5),LETTERS,letters),22,replace=TRUE),
collapse=''),
sep='') )
} else { #Call AFNI program
return(system('3dnewid -fun', intern=TRUE) );
}
}
AFNI.view2viewtype <- function(view="+orig"){
if (is.null(view)) return(0)
if (length(grep("orig", view, ignore.case=TRUE))) return(0)
if (length(grep("acpc", view, ignore.case=TRUE))) return(1)
if (length(grep("tlrc", view, ignore.case=TRUE))) return(2)
return(0)
}
AFNI.viewtype2view <- function(view=0){
if (is.null(view)) return("+orig")
if (view == 0) return("+orig")
if (view == 1) return("+acpc")
if (view == 2) return("+tlrc")
return("+orig")
}
dset.view <- function(dset) {
dd <- dset.attr(dset$NI_head, "SCENE_DATA")
return( AFNI.viewtype2view(dd[1]))
}
dset.gridmatch <- function(mdset, idset) {
if (is.null(mdset) || is.null(idset)) return(FALSE)
if (dset.nvox(mdset) != dset.nvox(idset)) return(FALSE)
return(TRUE);
}
write.c.AFNI <- function( filename, dset=NULL, label=NULL, space=NULL,
note=NULL, origin=NULL, delta=NULL,
orient=NULL,
idcode=NULL, defhead=NULL,
verb = 1,
maskinf=0, scale = TRUE,
overwrite=FALSE, addFDR=0,
statsym=NULL, view=NULL,
com_hist=NULL, TR=NULL, type=NULL) {
an <- parse.AFNI.name(filename);
if (verb > 1) {
show.AFNI.name(an);
}
if (class(dset) == "AFNI_c_dataset") {
#Assume all is in dset already
if (is.null(dset$NI_head)) {
err.AFNI("Old dset?");
return(NULL);
}
} else if (class(dset) == "AFNI_R_dataset") {
err.AFNI("Should not be with R dataset here")
return(NULL);
} else {
#Assume user sent in an array, not a dset
dset <- array2dset(dset,meth='clib')
}
if (!is.loaded('R_THD_write_dset')) {
err.AFNI("Missing R_io.so.");
return(NULL);
}
if (is.null(view)) {
if (is.null(defhead)) {
view <- "+tlrc"
}else {
view <- dset.view(defhead)
}
}
#Revert to old default if view is still not set
if (is.null(view)) view <- "+tlrc"
#Setup the basics (might be repetitive, oh well)
dset$NI_head <- dset.attr(dset$NI_head, "dim", val = dset.dimBRKarray(dset))
#Make sure TYPESTRING and SCENE_DATA are coherent
dset$NI_head <- dset.attr( dset$NI_head, "TYPESTRING",
default = '3DIM_GEN_FUNC')
dset$NI_head <- dset.attr( dset$NI_head, "SCENE_DATA",
default = c(AFNI.view2viewtype(view),11,3, rep(0,5)))
taxnums <- NULL;
taxfloats <- NULL;
taxoff <- NULL;
if (!is.null(defhead)) {
taxnums <- dset.attr( defhead, "TAXIS_NUMS");
taxfloats <- dset.attr( defhead, "TAXIS_FLOATS");
taxoff <- dset.attr( defhead, "TAXIS_OFFSETS");
}
#The user options
if (is.null(idcode)) idcode <- newid.AFNI();
if (!is.null(idcode)) dset$NI_head <-
dset.attr(dset$NI_head, "IDCODE_STRING", val=idcode)
if (!is.null(TR)) {
if (is.null(taxnums)) {
#77002 for seconds
taxnums <- c(dset.nval(dset), 0, 77002, -999, -999, -999, -999, -999);
}
if (is.null(taxfloats)) {
taxfloats <- c(0, TR, 0, 0, 0, -999999.000, -999999.000, -999999.0000);
}
if (is.null(taxoff)) {
taxoff <- rep(0, dset.nslize(dset));
}
dset$NI_head <- dset.attr( dset$NI_head, "TAXIS_NUMS", val = taxnums);
dset$NI_head <- dset.attr( dset$NI_head, "TAXIS_FLOATS", val = taxfloats);
dset$NI_head <- dset.attr( dset$NI_head, "TAXIS_OFFSETS", val = taxoff);
dset$NI_head <- dset.attr(dset$NI_head, "TR", val = TR);
} else { #use defhead values, setup from above
if (!is.null(taxnums)) {
dset$NI_head <- dset.attr( dset$NI_head, "TAXIS_NUMS", val = taxnums);
}
if (!is.null(taxfloats)) {
dset$NI_head <- dset.attr( dset$NI_head, "TAXIS_FLOATS",
val = taxfloats);
}
if (!is.null(taxoff)) {
dset$NI_head <- dset.attr( dset$NI_head, "TAXIS_OFFSETS",
val = taxoff);
}
}
if (is.null(origin) && !is.null(defhead))
origin <- dset.attr(defhead,"ORIGIN")
if (!is.null(origin)) dset$NI_head <-
dset.attr(dset$NI_head, "ORIGIN", val=origin)
if (is.null(delta) && !is.null(defhead))
delta <- dset.attr(defhead,"DELTA")
if (!is.null(delta)) dset$NI_head <-
dset.attr(dset$NI_head, "DELTA", val=delta)
IJK_TO_DICOM_REAL <- NULL
if (!is.null(defhead))
IJK_TO_DICOM_REAL <- dset.attr(defhead$NI_head,"IJK_TO_DICOM_REAL")
if (!is.null(IJK_TO_DICOM_REAL)) dset$NI_head <-
dset.attr(dset$NI_head, "IJK_TO_DICOM_REAL", val=IJK_TO_DICOM_REAL)
if (is.null(orient) && !is.null(defhead))
orient <- dset.attr(defhead,"ORIENT_SPECIFIC")
if (!is.null(orient)) dset$NI_head <-
dset.attr(dset$NI_head, "ORIENT_SPECIFIC",
val = orcode.AFNI(orient))
if (is.null(label) && !is.null(defhead))
label <- dset.attr(defhead,"BRICK_LABS")
if (!is.null(label)) dset$NI_head <-
dset.attr(dset$NI_head, "BRICK_LABS",
val = label);
if (is.null(space) && !is.null(defhead))
space <- dset.attr(defhead,"TEMPLATE_SPACE")
if (!is.null(space)) dset$NI_head <-
dset.attr(dset$NI_head, "TEMPLATE_SPACE",
val = space);
if (is.null(type) && !is.null(defhead)) {
tps = dset.attr(defhead,"BRICK_TYPES");
if (!is.null(tps) && length(tps)>0) type <- tps[1]
}
if (!is.null(type)) {
if (!is.null(type <- typecode.AFNI(type))) {
type <- rep(type, dset.dimBRKarray(dset)[4])
dset$NI_head <- dset.attr(dset$NI_head, "BRICK_TYPES", val=type)
}
}
if (is.null(note) && !is.null(defhead))
note <- dset.attr(defhead,"note")
if (!is.null(note)) dset$NI_head <-
dset.attr(dset$NI_head, "note",
val = note);
if (is.null(statsym) && !is.null(defhead))
statsym <- dset.attr(defhead,"statsym")
if (!is.null(statsym)) dset$NI_head <-
dset.attr(dset$NI_head, "statsym",
val = statsym);
if (maskinf) {
if (verb) note.AFNI("Masking infs");
dset$brk[!is.finite(dset$brk)]=0
}
uopts = list( scale = scale, overwrite=overwrite, debug=verb,
addFDR=addFDR, hist=com_hist);
dset$head <- unparse.c.AFNI.head(dset$NI_head)
rso <- .Call("R_THD_write_dset", name = as.character(filename),
dset = as.list(dset), opts = uopts)
invisible(dset);
}
write.AFNI <- function( filename, brk=NULL, label=NULL, space=NULL,
note=NULL, origin=NULL, delta=NULL,
orient=NULL,
idcode=NULL, defhead=NULL,
verb = 0,
maskinf=0, scale = TRUE,
meth='AUTO', addFDR=FALSE,
statsym=NULL, view=NULL,
com_hist=NULL, type=NULL, TR = NULL,
overwrite=FALSE) {
if (meth == 'AUTO') {
if (class(brk) == 'AFNI_R_dataset') {
meth = 'Rlib'
} else if (class(brk) == 'AFNI_c_dataset') {
meth = 'clib'
} else if (have_R_io()) {
meth <- 'clib'
} else {
meth <- 'Rlib'
}
}
an <- parse.AFNI.name(filename);
if (verb>1) {
show.AFNI.name(an);
}
if (an$type == "1D" || an$type == "Rs" || an$type == "1Ds") {
return(write.AFNI.matrix(drop(brk), an$orig_name))
}
if (meth == 'clib' || an$type == 'NIML' || class(brk) == "AFNI_c_dataset") {
return(write.c.AFNI(filename, dset=brk, label=label, space=space,
note=note, origin=origin, delta=delta,
orient=orient,
idcode=idcode, defhead=defhead,
verb = verb,
maskinf=maskinf, scale = scale, addFDR=addFDR,
statsym=statsym, view=view, com_hist=com_hist,
type=type, TR=TR, overwrite=overwrite))
}
#Switch to old default for view if not using write.c.AFNI
if (is.null(view)) view <- "+tlrc"
if (is.na(an$view)) {
err.AFNI('Bad filename for old writing method. Need view');
show.AFNI.name(an);
return(0)
}
if (addFDR && verb) {
warn.AFNI("addFDR is not supported for AFNI_R_dataset structs");
}
if (!is.null(statsym) && verb) {
warn.AFNI("statsym is not supported for AFNI_R_dataset structs");
}
#Call with brk if dset is sent in
if (class(brk) == "AFNI_R_dataset") {
return(write.AFNI( filename, brk$brk, label, space, note, origin, delta,
orient, idcode, defhead, verb, maskinf, scale))
}
#Set the defaults.
if (is.null(note)) note <- '';
if (is.null(idcode)) idcode <- newid.AFNI();
if (is.null(brk)) {
err.AFNI("NULL input");
return(0)
}
#make sure we don't have 3 dims array, need always 4d
ddb <- dimBRKarray(brk)
dim(brk) <- ddb
if (is.null(defhead)) { # No default header
if (verb) note.AFNI("Have no default header");
if (is.null(label)) label <- paste(c(1:ddb[4]),collapse='~');
if (is.null(space)) space <-'TLRC';
if (is.null(origin)) origin <- c(0,0,0)
if (is.null(delta)) delta <- c(4,4,4)
if (is.null(orient)) orient <- 'RAI'
} else { #When possible, call on default header
if (verb) note.AFNI("Have default header");
if (is.null(label)) {
if (!is.null(defhead$BRICK_LABS)) {
label <- gsub("^'", '', defhead$BRICK_LABS);
} else {
label <- paste(c(1:ddb[4]),collapse='~');
}
}
if (is.null(space)) {
if (!is.null(defhead$NI_head$TEMPLATE_SPACE$dat)) {
space <- defheadTEMPLATE_SPACE$dat;
} else {
space <- 'TLRC';
}
}
if (is.null(origin)) {
if (!is.null(defhead$ORIGIN)) {
origin <- defhead$ORIGIN;
} else {
origin <- c(0,0,0);
}
}
if (is.null(delta)) {
if (!is.null(defhead$DELTA)) {
delta <- defhead$DELTA;
} else {
delta <- c(4,4,4);
}
}
if (is.null(orient)) {
if (!is.null(defhead$ORIENT_SPECIFIC)) {
orient <- defhead$ORIENT_SPECIFIC;
} else {
orient <- 'RAI';
}
}
}
if (is.null(brk) || !is.numeric(brk)) {
err.AFNI("data array improperly formatted");
return(0);
}
if (verb) {
note.AFNI("Writing header")
}
conhead <- file(head.AFNI.name(an), "w")
writeChar(AFNIheaderpart("string-attribute","TEMPLATE_SPACE",space),
conhead,eos=NULL)
writeChar(AFNIheaderpart("string-attribute","HISTORY_NOTE",note),
conhead,eos=NULL)
writeChar(AFNIheaderpart("string-attribute","TYPESTRING","3DIM_HEAD_FUNC"),
conhead,eos=NULL)
writeChar(AFNIheaderpart("string-attribute","IDCODE_STRING",''),
conhead,eos=NULL)
writeChar(AFNIheaderpart("string-attribute","IDCODE_DATE",date()),
conhead,eos=NULL)
if (an$view == '+orig') { sv <- 0 }
else if (an$view == '+acpc') { sv <- 1 }
else if (an$view == '+tlrc') { sv <- 2 }
else { sv <- 0 }
writeChar(AFNIheaderpart("integer-attribute","SCENE_DATA",
c(sv,11,1,-999,-999,-999,-999,-999)),
conhead,eos=NULL)
writeChar(AFNIheaderpart("integer-attribute","ORIENT_SPECIFIC",
orcode.AFNI(orient)),
conhead,eos=NULL)
writeChar(AFNIheaderpart("float-attribute","ORIGIN",origin),
conhead,eos=NULL)
writeChar(AFNIheaderpart("float-attribute","DELTA",delta),
conhead,eos=NULL)
if (maskinf) {
if (verb) note.AFNI("Masking infs");
brk[!is.finite(brk)]=0
}
mm <- minmax(brk)
if (verb) {
note.AFNI("minmax");
str(mm)
}
writeChar(AFNIheaderpart("float-attribute","BRICK_STATS",mm),
conhead,eos=NULL)
writeChar(AFNIheaderpart("integer-attribute","DATASET_RANK",
c(3,ddb[4],0,0,0,0,0,0)),
conhead,eos=NULL)
writeChar(AFNIheaderpart("integer-attribute","DATASET_DIMENSIONS",
c(ddb[1:3],0,0)),
conhead,eos=NULL)
if (is.null(type)) type <- 1
writeChar(AFNIheaderpart("integer-attribute","BRICK_TYPES",
rep(typecode.AFNI(type),ddb[4])),
conhead,eos=NULL)
if (scale) {
if (verb) {
note.AFNI("Computing scaling factors");
}
scale_fac <- rep(0,ddb[4])
for (k in 1:ddb[4]) {
scale_fac[k] <- max(abs(mm[2*k-1]),abs(mm[2*k]))/32767
}
} else {
scale_fac <- rep(0,ddb[4])
}
if (verb) {
note.AFNI("scaling factors:");
str(scale_fac);
}
writeChar(AFNIheaderpart("float-attribute","BRICK_FLOAT_FACS",scale_fac),
conhead,eos=NULL)
writeChar(AFNIheaderpart("string-attribute","BRICK_LABS",
paste(label,collapse="~")),
conhead,eos=NULL)
writeChar(AFNIheaderpart("string-attribute","BYTEORDER_STRING","MSB_FIRST"),
conhead,eos=NULL)
#The tross logger
trlog <- Sys.getenv("AFNI_HISTDB_SCRIPT")
if ( trlog != '') {
writeChar(AFNIheaderpart("string-attribute","HISTDB_SCRIPT", trlog ),
conhead,eos=NULL) ;
}
close(conhead)
# Write BRIK
conbrik <- file(brik.AFNI.name(an), "wb")
if (0) { #ZSS: old method,
# runs out of memory for large dsets
# when scaling. Code kept here for testing
if (scale) {
if (verb) {
note.AFNI("Applying scaling ")
}
for (k in 1:ddb[4]) {
brk[,,,k] <- brk[,,,k] / scale_fac[k]
}
}
dim(brk) <- NULL #Don't know why this was done here ....
if (verb) {
note.AFNI("Writing brik")
}
writeBin(as.integer(brk), conbrik,size=2, endian="big")
} else {
if (verb) {
note.AFNI("Writing /Scaling, meth 2")
}
#Write on sub-brick at a time to reduce memory use.
# as.integer will allocate a new copy
if (scale) {
for (k in 1:ddb[4]) {
writeBin(as.integer(brk[,,,k] / scale_fac[k]),
conbrik,size=2, endian="big")
}
} else {
for (k in 1:ddb[4]) {
writeBin(as.integer(brk[,,,k]),
conbrik,size=2, endian="big")
}
}
}
close(conbrik)
return(1);
}
read.NIFTI <- function(filename) {
fileparts <- strsplit(filename,"\\.")[[1]]
ext <- tolower(fileparts[length(fileparts)])
if (ext == "nii") {
filename.nii <- filename
filename.hdr <- paste(c(fileparts[-length(fileparts)],"hdr"),collapse=".")
filename.img <- paste(c(fileparts[-length(fileparts)],"img"),collapse=".")
} else if (ext == "hdr") {
filename.hdr <- filename
filename.img <- paste(c(fileparts[-length(fileparts)],"img"),collapse=".")
} else if (ext == "img") {
filename.hdr <- paste(c(fileparts[-length(fileparts)],"hdr"),collapse=".")
filename.img <- filename
} else {
filename.nii <- paste(filename,".nii",sep="")
filename.hdr <- paste(filename,".hdr",sep="")
filename.img <- paste(filename,".img",sep="")
}
if ((ext != "hdr") && (ext != "img") && (!is.na(file.info(filename.nii)$size))) {
con <- file(filename.nii,"rb")
header <- read.NIFTI.header(con)
if (!(header$magic == "n+1") && !(header$magic == "ni1"))
warning("Hmmm! Dont see the magic NIFTI string! Try to proceed, but maybe some weird results will occur!");
bytes <- header$voxoffset - 348
header$extension <- readBin(con,"raw",bytes)
} else {
if (is.na(file.info(filename.hdr)$size) | (file.info(filename.hdr)$size < 348))
stop("Hmmm! This does not seem to be a NIFTI header (hdr/img-pair)! Wrong size or does not exist!");
con <- file(filename.hdr,"rb")
header <- read.NIFTI.header(con)
header$extension <- NULL
close(con)
if (is.na(file.info(filename.img)$size))
stop("Hmmm! This does not seem to be a NIFTI header (hdr/img-pair)! img-file not found!");
con <- file(filename.img,"rb")
}
dx <- header$dimension[2]
dy <- header$dimension[3]
dz <- header$dimension[4]
dt <- header$dimension[5]
endian <- header$endian
if (header$datatype == 1) { # logical
what <- "raw"
signed <- TRUE
size <- 1
} else if (header$datatype == 2) { # unsigned char????
what <- "int"
signed <- FALSE
size <- if (header$bitpix) header$bitpix/8 else 2
} else if (header$datatype == 4) { # signed short
what <- "int"
signed <- TRUE
size <- if (header$bitpix) header$bitpix/8 else 2
} else if (header$datatype == 8) { # signed integer
what <- "int"
signed <- TRUE
size <- if (header$bitpix) header$bitpix/8 else 4
} else if (header$datatype == 16) { # float
what <- "double"
signed <- TRUE
size <- if (header$bitpix) header$bitpix/8 else 4
} else if (header$datatype == 32) { # complex
what <- "complex"
signed <- TRUE
size <- if (header$bitpix) header$bitpix/8 else 8
} else if (header$datatype == 64) { # double
what <- "double"
signed <- TRUE
size <- if (header$bitpix) header$bitpix/8 else 8
} else { # all other
what <- "raw"
signed <- TRUE
size <- 1
}
brk <- readBin(con, what, n=dx*dy*dz*dt, size=size, signed=signed, endian=endian)
close(con)
if (min(abs(header$pixdim[2:4])) != 0) {
weights <-
abs(header$pixdim[2:4]/min(abs(header$pixdim[2:4])))
} else {
weights <- NULL
}
dim(brk) <- c(dx,dy,dz,dt)
mask <- array(TRUE,c(dx,dy,dz))
mask[brk[,,,1] < quantile(brk[,,,1],0.75)] <- FALSE
z <- list(brk=writeBin(as.numeric(brk),raw(),4),format="NIFTI",delta=header$pixdim[2:4],
origin=NULL,orient=NULL,dim=header$dimension[2:5],weights=weights,header=header,mask=mask)
class(z) <- "AFNI_R_dataset"
invisible(z)
}
extract.data <- function(z,what="data") {
if (!("AFNI_R_dataset"%in%class(z)) || !("AFNI_c_dataset"%in%class(z))) {
warning("extract.data: data not of class <AFNI_[cR]_dataset>. Try to proceed but strange things may happen")
}
if (what=="residuals") {
if(!is.null(z$resscale)){
brk <- readBin(z$res,"integer",prod(z$dim),2)*z$resscale
dim(brk) <- z$dim
} else {
warning("extract.data: No residuals available, returning NULL")
brk <- NULL
}
} else {
brk <- readBin(z$brk,"numeric",prod(z$dim),4)
dim(brk) <- z$dim
}
invisible(brk)
}
# GC 01/12/2015: Used in 3dMVM.T and 3dLME.R to construct GLT and GLF stuff
#gl_Constr2 <- function(n_gl, code, lop) { # n_gl: number of tests: lop$num_glt or lop$num_glf; code: lop$glfCode
# if (n_gl > 0) {
# outList <- vector('list', 3)
# outList[[1]] <- vector('list', n_gl)
# outList[[2]] <- vector('list', n_gl)
# outList[[3]] <- vector('list', n_gl)
# for (n in 1:n_gl) { # assuming each GLT has one slope involved and placed last
# # if(length(lop$QV)==0) outList[[1]][[n]] <- glfConstr(code[[n]], lop$dataStr) else {
# if((length(lop$QV)==0) & is.na(lop$vVars)) outList[[1]][[n]] <- glfConstr(code[[n]], lop$dataStr) else {
# if((length(lop$QV)>0) & any(lop$QV %in% code[[n]])) {
# QVpos <- which(code[[n]] %in% lop$QV)
# if(is.na(code[[n]][QVpos+2])) { # test for covariate effect
# if(QVpos==1) outList[[1]][[n]] <- NA else
# outList[[1]][[n]] <-glfConstr(code[[n]][-c(QVpos, QVpos+1)], lop$dataStr)
# outList[[2]][[n]] <- code[[n]][QVpos]
# } else { # effect at a specific covariate value
# if(QVpos==1) outList[[1]][[n]] <- NA else
# outList[[1]][[n]] <-glfConstr(code[[n]][-(QVpos:(QVpos+2))], lop$dataStr)
# outList[[3]][[n]] <- as.numeric(code[[n]][QVpos+2])
# names(outList[[3]][[n]]) <- code[[n]][QVpos]
# } # if(is.na(lop$gltCode[[n]][QVpos+2]))
# } else if(!is.na(lop$vVars) & any(lop$vQV %in% code[[n]])) { # voxel-wise covariate
# vQVpos <- which(code[[n]] %in% lop$vQV)
# if(is.na(code[[n]][vQVpos+2])) { # test for covariate effect
# if(vQVpos==1) outList[[1]][[n]] <- NA else
# outList[[1]][[n]] <-glfConstr(code[[n]][-c(vQVpos, vQVpos+1)], lop$dataStr)
# outList[[2]][[n]] <- code[[n]][vQVpos]
# } else { # effect at a specific covariate value
# if(vQVpos==1) outList[[1]][[n]] <- NA else
# outList[[1]][[n]] <-glfConstr(code[[n]][-(vQVpos:(vQVpos+2))], lop$dataStr)
# outList[[3]][[n]] <- as.numeric(code[[n]][vQVpos+2])
# names(outList[[3]][[n]]) <- code[[n]][vQVpos]
# } # if(is.na(lop$gltCode[[n]][vQVpos+2]))
# } else outList[[1]][[n]] <- glfConstr(code[[n]], lop$dataStr) # if((length(lop$QV)>0) & any(lop$QV %in% lop$gltCode[[n]]))
# }
# }
# }
# return(outList)
#}
gl_Constr <- function(n_gl, code, lop) { # n_gl: number of tests: lop$num_glt or lop$num_glf; code: lop$glfCode
if (n_gl > 0) {
outList <- vector('list', 3)
outList[[1]] <- vector('list', n_gl)
outList[[2]] <- vector('list', n_gl)
outList[[3]] <- vector('list', n_gl)
comQV <- c(lop$QV, lop$vVars)
for (n in 1:n_gl) { # assuming each GLT has one slope involved and placed last
if(length(comQV)==0) outList[[1]][[n]] <- glfConstr(code[[n]], lop$dataStr) else {
if((length(comQV)>0) & any(comQV %in% code[[n]])) {
QVpos <- which(code[[n]] %in% comQV)
if(is.na(code[[n]][QVpos+2])) { # test for covariate effect
if(QVpos==1) outList[[1]][[n]] <- NA else
outList[[1]][[n]] <-glfConstr(code[[n]][-c(QVpos, QVpos+1)], lop$dataStr)
outList[[2]][[n]] <- code[[n]][QVpos]
} else { # effect at a specific covariate value
if(QVpos==1) outList[[1]][[n]] <- NA else
outList[[1]][[n]] <-glfConstr(code[[n]][-(QVpos:(QVpos+2))], lop$dataStr)
outList[[3]][[n]] <- as.numeric(code[[n]][QVpos+2])
names(outList[[3]][[n]]) <- code[[n]][QVpos]
} # if(is.na(lop$gltCode[[n]][QVpos+2]))
} else outList[[1]][[n]] <- glfConstr(code[[n]], lop$dataStr) # if((length(lop$QV)>0) & any(lop$QV %in% lop$gltCode[[n]]))
}
if(is.null(outList[[1]][[n]])) errex.AFNI(paste("Inappropriate coding in test No.", n, "! \n ", sep = ""))
}
}
return(outList)
}
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.