4
# Copyright (c) 2000-2007 Gary Strangman. All rights reserved
8
# This software is provided "as-is". There are no expressed or implied
9
# warranties of any kind, including, but not limited to, the warranties
10
# of merchantability and fittness for a given application. In no event
11
# shall Gary Strangman be liable for any direct, indirect, incidental,
12
# special, exemplary or consequential damages (including, but not limited
13
# to, loss of use, data or profits, or business interruption) however
14
# caused and on any theory of liability, whether in contract, strict
15
# liability or tort (including negligence or otherwise) arising in any way
16
# out of the use of this software, even if advised of the possibility of
19
# Comments and/or additions are welcome (send e-mail to:
20
# strang@nmr.mgh.harvard.edu).
23
Defines a number of functions for pseudo-command-line OS functionality.
26
pwd <-- can be used WITHOUT parens
29
get(namepatterns,verbose=1)
30
getstrings(namepatterns,verbose=1)
31
put(outlist,filename,writetype='w')
32
aget(namepatterns,verbose=1)
33
aput(outarray,filename,writetype='w')
34
bget(filename,numslices=1,xsize=64,ysize=64)
36
bput(outarray,filename,writeheader=0,packstring='h',writetype='wb')
43
## 07-11-26 ... more numpy conversion work
44
## 06-08-07 ... converted to numpy, changed version to 0.6
45
## 06-02-03 ... added add2afnihistory() to load modify afni HISTORY_NOTEs,
46
## and added that option to array2afni output
47
## 04-06-14 ... added getafniparam() to load in afni values from HEAD files
48
## 03-04-09 ... fixed brikget to load based on datatype, changed array2afni
49
## so that the default sliceorder is altplus, not seqplus
50
## 02-11-20 ... added binget(), binput(), array2afni(), version 0.5
51
## 02-10-20 ... added find_dirs() function, changed version to 0.4
52
## 01-11-15 ... changed aput() and put() to accept a delimiter
53
## 01-04-19 ... added oneperline option to put() function
54
## 99-11-07 ... added DAs quick flat-text-file loaders, load() and fload()
55
## 99-11-01 ... added version number (0.1) for distribution
56
## 99-08-30 ... Put quickload in here
57
## 99-06-27 ... Changed bget thing back ... confused ...
58
## 99-06-24 ... exchanged xsize and ysize in bget for non-square images (NT??)
59
## modified bget to raise an IOError when file not found
60
## 99-06-12 ... added load() and save() aliases for aget() and aput() (resp.)
61
## 99-04-13 ... changed aget() to ignore (!!!!) lines beginning with # or %
62
## 99-01-17 ... changed get() so ints come in as ints (not floats)
73
import glob, re, string, types, os, struct, copy, time, tempfile, sys
81
Wraps a function so that if it's entered *by itself*
82
in the interpreter without ()'s, it gets called anyway
85
def __init__(self, f):
98
Changes the working python directory for the interpreter.
100
Usage: cd(directory) where 'directory' is a string
108
Changes the working python directory for the interpreter.
110
Usage: pwd (no parens needed)
118
Produces a directory listing. Default is the current directory.
126
def rename(source, dest):
128
Renames files specified by UNIX inpattern to those specified by UNIX
129
outpattern. Can only handle a single '*' in the two patterns!!!
131
Usage: rename(source, dest) e.g., rename('*.txt', '*.c')
133
infiles = glob.glob(source)
135
incutindex = string.index(source,'*')
136
outcutindex = string.index(source,'*')
137
findpattern1 = source[0:incutindex]
138
findpattern2 = source[incutindex+1:]
139
replpattern1 = dest[0:incutindex]
140
replpattern2 = dest[incutindex+1:]
141
for fname in infiles:
143
newname = re.sub(findpattern1,replpattern1,fname,1)
144
if outcutindex < len(dest)-1:
146
lastone = string.rfind(newname,replpattern2)
147
newname = newname[0:lastone] + re.sub(findpattern2,replpattern2,fname[lastone:],1)
149
lastone = string.rfind(fname,findpattern2)
151
newname = fname[0:lastone]
152
newname = newname + re.sub(findpattern2,replpattern2,fname[lastone:],1)
154
os.rename(fname,newname)
158
def get(namepatterns,verbose=1):
160
Loads a list of lists from text files (specified by a UNIX-style
161
wildcard filename pattern) and converts all numeric values to floats.
162
Uses the glob module for filename pattern conversion. Loaded filename
163
is printed if verbose=1.
165
Usage: get(namepatterns,verbose=1)
166
Returns: a 1D or 2D list of lists from whitespace delimited text files
167
specified by namepatterns; numbers that can be converted to floats
171
if type(namepatterns) in [ListType,TupleType]:
172
for item in namepatterns:
173
fnames = fnames + glob.glob(item)
175
fnames = glob.glob(namepatterns)
179
print 'NO FILENAMES MATCH ('+namepatterns+') !!'
183
print fnames # so user knows what has been loaded
185
for i in range(len(fnames)):
186
file = open(fnames[i])
187
newelements = map(string.split,file.readlines())
188
for i in range(len(newelements)):
189
for j in range(len(newelements[i])):
191
newelements[i][j] = string.atoi(newelements[i][j])
194
newelements[i][j] = string.atof(newelements[i][j])
197
elements = elements + newelements
198
if len(elements)==1: elements = elements[0]
202
def getstrings(namepattern,verbose=1):
204
Loads a (set of) text file(s), with all elements left as string type.
205
Uses UNIX-style wildcards (i.e., function uses glob). Loaded filename
206
is printed if verbose=1.
208
Usage: getstrings(namepattern, verbose=1)
209
Returns: a list of strings, one per line in each text file specified by
212
fnames = glob.glob(namepattern)
215
print 'NO FILENAMES MATCH ('+namepattern+') !!'
220
for filename in fnames:
221
file = open(filename)
222
newelements = map(string.split,file.readlines())
223
elements = elements + newelements
227
def put(outlist,fname,writetype='w',oneperline=0,delimit=' '):
229
Writes a passed mixed-type list (str and/or numbers) to an output
230
file, and then closes the file. Default is overwrite the destination
233
Usage: put(outlist,fname,writetype='w',oneperline=0,delimit=' ')
236
if type(outlist) in [N.ndarray]:
237
aput(outlist,fname,writetype)
239
if type(outlist[0]) not in [ListType,TupleType]: # 1D list
240
outfile = open(fname,writetype)
242
outlist = pstat.list2string(outlist,delimit)
243
outfile.write(outlist)
245
else: # they want one element from the list on each file line
247
outfile.write(str(item)+'\n')
249
else: # 2D list (list-of-lists)
250
outfile = open(fname,writetype)
252
outfile.write(pstat.list2string(row,delimit))
259
if type(x)==StringType:
266
def aget(namepattern,verbose=1):
268
Loads an array from 2D text files (specified by a UNIX-style wildcard
269
filename pattern). ONLY 'GET' FILES WITH EQUAL NUMBERS OF COLUMNS
270
ON EVERY ROW (otherwise returned array will be zero-dimensional).
272
Usage: aget(namepattern)
273
Returns: an array of integers, floats or objects (type='O'), depending on the
274
contents of the files specified by namepattern
276
fnames = glob.glob(namepattern)
279
print 'NO FILENAMES MATCH ('+namepattern+') !!'
284
for filename in fnames:
285
file = open(filename)
286
newelements = file.readlines()
288
for row in range(len(newelements)):
289
if (newelements[row][0]=='%' or newelements[row][0]=='#'
290
or len(newelements[row])==1 or newelements[row][0]=='\r'):
295
newelements = map(string.split,newelements)
296
for i in range(len(newelements)):
297
for j in range(len(newelements[i])):
299
newelements[i][j] = string.atof(newelements[i][j])
302
elements = elements + newelements
303
for row in range(len(elements)):
304
if N.add.reduce(N.array(map(isstring,elements[row])))==len(elements[row]):
305
print "A row of strings was found. Returning a LIST."
308
elements = N.array(elements)
310
elements = N.array(elements,dtype='O')
314
def aput(outarray,fname,writetype='w',delimit=' '):
316
Sends passed 1D or 2D array to an output file and closes the file.
318
Usage: aput(outarray,fname,writetype='w',delimit=' ')
321
outfile = open(fname,writetype)
322
if len(outarray.shape) == 1:
323
outarray = outarray[N.newaxis,:]
324
if len(outarray.shape) > 2:
325
raise TypeError, "put() and aput() require 1D or 2D arrays. Otherwise use some kind of pickling."
326
else: # must be a 2D array
328
outfile.write(string.join(map(str,row),delimit))
334
def bget(imfile,shp=None,unpackstr=N.int16,bytesperpixel=2.0,sliceinit=0):
336
Reads in a binary file, typically with a .bshort or .bfloat extension.
337
If so, the last 3 parameters are set appropriately. If not, the last 3
338
parameters default to reading .bshort files (2-byte integers in big-endian
341
Usage: bget(imfile,shp=None,unpackstr=N.int16,bytesperpixel=2.0,sliceinit=0)
343
if imfile[:3] == 'COR':
344
return CORget(imfile)
345
if imfile[-2:] == 'MR':
346
return mrget(imfile,unpackstr)
347
if imfile[-4:] == 'BRIK':
348
return brikget(imfile,unpackstr,shp)
349
if imfile[-3:] in ['mnc','MNC','inc','INC']:
350
return mincget(imfile,unpackstr,shp)
351
if imfile[-3:] == 'img':
352
return mghbget(imfile,unpackstr,shp)
353
if imfile[-6:] == 'bshort' or imfile[-6:] == 'bfloat':
355
return mghbget(imfile,unpackstr=unpackstr,bytesperpixel=bytesperpixel,sliceinit=sliceinit)
357
return mghbget(imfile,shp[0],shp[1],shp[2],unpackstr,bytesperpixel,sliceinit)
362
Reads a binary COR-nnn file (flattening file).
364
Usage: CORget(imfile)
365
Returns: 2D array of 16-bit ints
367
d=braw(infile,N.int8)
369
d = N.where(d>=0,d,256+d)
373
def mincget(imfile,unpackstr=N.int16,shp=None):
375
Loads in a .MNC file.
377
Usage: mincget(imfile,unpackstr=N.int16,shp=None) default shp = -1,20,64,64
381
os.system('mincextract -short -range 0 4095 -image_range 0 4095 ' +
382
imfile+' > minctemp.bshort')
384
d = braw('minctemp.bshort',unpackstr)
386
print "Couldn't find file: "+imfile
387
raise IOError, "Couldn't find file in mincget()"
391
os.system('rm minctemp.bshort')
395
def brikget(imfile,unpackstr=N.int16,shp=None):
397
Gets an AFNI BRIK file.
399
Usage: brikget(imfile,unpackstr=N.int16,shp=None) default shp: (-1,48,61,51)
404
file = open(imfile, "rb")
406
print "Couldn't find file: "+imfile
407
raise IOError, "Couldn't find file in brikget()"
409
header = imfile[0:-4]+'HEAD'
410
lines = open(header).readlines()
411
for i in range(len(lines)):
412
if string.find(lines[i],'DATASET_DIMENSIONS') <> -1:
413
dims = string.split(lines[i+2][0:string.find(lines[i+2],' 0')])
414
dims = map(string.atoi,dims)
415
if string.find(lines[i],'BRICK_FLOAT_FACS') <> -1:
416
count = string.atoi(string.split(lines[i+1])[2])
418
for j in range(int(N.ceil(count/5.))):
419
mults += map(string.atof,string.split(lines[i+2+j]))
420
mults = N.array(mults)
421
if string.find(lines[i],'BRICK_TYPES') <> -1:
423
first5 = map(string.atoi,string.split(first5))
429
unpackstr = N.float32
431
unpackstr = N.complex32
435
print "No header file. Continuing ..."
439
print 'Using unpackstr:',unpackstr #,', bytesperpixel=',bytesperpixel
441
file = open(imfile, "rb")
444
# the > forces big-endian (for or from Sun/SGI)
445
bdata = N.fromstring(bdata,unpackstr)
446
# littleEndian = ( struct.pack('i',1)==struct.pack('<i',1) )
447
if (max(bdata)>1e30):
448
bdata = bdata.byteswap()
452
print 'Incorrect shape ...',shp,len(bdata)
453
raise ValueError, 'Incorrect shape for file size'
457
if N.sum(mults) == 0:
460
multshape = [1]*len(bdata.shape)
461
for i in range(len(bdata.shape)):
462
if len(mults) == bdata.shape[i]:
463
multshape[i] = len(mults)
465
mults.shape = multshape
470
def mghbget(imfile,numslices=-1,xsize=64,ysize=64,
471
unpackstr=N.int16,bytesperpixel=2.0,sliceinit=0):
473
Reads in a binary file, typically with a .bshort or .bfloat extension.
474
If so, the last 3 parameters are set appropriately. If not, the last 3
475
parameters default to reading .bshort files (2-byte integers in big-endian
478
Usage: mghbget(imfile, numslices=-1, xsize=64, ysize=64,
479
unpackstr=N.int16, bytesperpixel=2.0, sliceinit=0)
482
file = open(imfile, "rb")
484
print "Couldn't find file: "+imfile
485
raise IOError, "Couldn't find file in bget()"
487
header = imfile[0:-6]+'hdr'
488
vals = get(header,0) # '0' means no missing-file warning msg
489
if type(vals[0]) == ListType: # it's an extended header
490
xsize = int(vals[0][0])
491
ysize = int(vals[0][1])
492
numslices = int(vals[0][2])
496
numslices = int(vals[2])
498
print "No header file. Continuing ..."
501
if suffix == 'bshort':
503
elif suffix[-3:] == 'img':
505
elif suffix == 'bfloat':
506
unpackstr = N.float32
510
print 'Not a bshort, bfloat or img file.'
511
print 'Using unpackstr:',unpackstr,', bytesperpixel=',bytesperpixel
514
file = open(imfile, "rb")
517
numpixels = len(bdata) / bytesperpixel
519
raise ValueError, "Incorrect file size in fmri.bget()"
520
else: # the > forces big-endian (for or from Sun/SGI)
521
bdata = N.fromstring(bdata,unpackstr)
522
# littleEndian = ( struct.pack('i',1)==struct.pack('<i',1) )
524
# bdata = bdata.byteswap()
525
if (max(bdata)>1e30):
526
bdata = bdata.byteswap()
527
if suffix[-3:] == 'img':
529
numslices = len(bdata)/8200 # 8200=(64*64*2)+8 bytes per image
532
slices = N.zeros((numslices,xsize,ysize),N.int32)
533
for i in range(numslices):
534
istart = i*8 + i*xsize*ysize
535
iend = i*8 + (i+1)*xsize*ysize
537
slices[i] = N.reshape(N.array(bdata[istart:iend]),(xsize,ysize))
540
slices = N.reshape(N.array(bdata),[xsize,ysize])
542
slices = N.reshape(N.array(bdata),[numslices,xsize,ysize])
548
def braw(fname,btype,shp=None):
550
Opens a binary file, unpacks it, and returns a flat array of the
551
type specified. Use Numeric types ... N.float32, N.int64, etc.
553
Usage: braw(fname,btype,shp=None)
554
Returns: flat array of floats, or ints (if btype=N.int16)
556
file = open(fname,'rb')
558
bdata = N.fromstring(bdata,btype)
559
# littleEndian = ( struct.pack('i',1)==struct.pack('<i',1) )
561
# bdata = bdata.byteswap() # didn't used to need this with '>' above
562
if (max(bdata)>1e30):
563
bdata = bdata.byteswap()
570
return N.array(bdata)
573
def glget(fname,btype):
575
Load in a file containing pixels from glReadPixels dump.
577
Usage: glget(fname,btype)
578
Returns: array of 'btype elements with shape 'shape', suitable for im.ashow()
580
d = braw(fname,btype)
585
shp = N.fromstring(shp,N.int32)
586
shp[0],shp[1] = shp[1],shp[0]
588
carray = N.reshape(d,shp)
599
carray = N.array([r,g,b])
601
outstr = "glget: shape not correct for data of length "+str(len(d))
602
raise ValueError, outstr
606
def mget(fname,btype):
608
Load in a file that was saved from matlab
610
Usage: mget(fname,btype)
612
d = braw(fname,btype)
614
header = fname[0:-6]+'hdr'
615
vals = get(header,0) # '0' means no missing-file warning msg
616
if type(vals[0]) == ListType: # it's an extended header
617
xsize = int(vals[0][0])
618
ysize = int(vals[0][1])
619
numslices = int(vals[0][2])
623
numslices = int(vals[2])
624
print xsize,ysize,numslices, d.shape
626
print "No header file. Continuing ..."
628
d.shape = [ysize,xsize]
629
return N.transpose(d)*1
631
d.shape = [numslices,ysize,xsize]
632
return N.transpose(d)*1
635
def mput(outarray,fname,writeheader=0,btype=N.int16):
637
Save a file for use in matlab.
639
outarray = N.transpose(outarray)
640
outdata = N.ravel(outarray).astype(btype)
641
outdata = outdata.tostring()
642
outfile = open(fname,'wb')
643
outfile.write(outdata)
647
suffixindex = string.rfind(fname,'.')
648
hdrname = fname[0:suffixindex]
651
if len(outarray.shape) == 2:
652
hdr = [outarray.shape[1],outarray.shape[0], 1, 0]
654
hdr = [outarray.shape[2],outarray.shape[1],outarray.shape[0], 0,'\n']
656
outfile = open(hdrname+'.hdr','w')
657
outfile.write(pstat.list2string(hdr))
662
def bput(outarray,fname,writeheader=0,packtype=N.int16,writetype='wb'):
664
Writes the passed array to a binary output file, and then closes
665
the file. Default is overwrite the destination file.
667
Usage: bput (outarray,filename,writeheader=0,packtype=N.int16,writetype='wb')
670
if suffix == 'bshort':
672
elif suffix == 'bfloat':
675
print 'Not a bshort or bfloat file. Using packtype=',packtype
677
outdata = N.ravel(outarray).astype(packtype)
678
# littleEndian = ( struct.pack('i',1)==struct.pack('<i',1) )
680
# outdata = outdata.byteswap()
681
outdata = outdata.tostring()
682
outfile = open(fname,writetype)
683
outfile.write(outdata)
687
suffixindex = string.rfind(fname,'.')
688
hdrname = fname[0:suffixindex]
691
if len(outarray.shape) == 2:
692
hdr = [outarray.shape[0],outarray.shape[1], 1, 0]
694
hdr = [outarray.shape[1],outarray.shape[2],outarray.shape[0], 0,'\n']
696
outfile = open(hdrname+'.hdr','w')
697
outfile.write(pstat.list2string(hdr))
702
def mrget(fname,datatype=N.int16):
704
Opens a binary .MR file and clips off the tail data portion of it, returning
705
the result as an array.
707
Usage: mrget(fname,datatype=N.int16)
709
d = braw(fname,datatype)
711
return N.reshape(d[-512*512:],(512,512))
712
elif len(d) > 320*320:
713
return N.reshape(d[-320*320:],(320,320))
714
elif len(d) > 256*256:
715
return N.reshape(d[-256*256:],(256,256))
716
elif len(d) > 128*128:
717
return N.reshape(d[-128*128:],(128,128))
719
return N.reshape(d[-64*64:],(64,64))
721
return N.reshape(d[-32*32:],(32,32))
724
def quickload(fname,linestocut=4):
726
Quickly loads in a long text file, chopping off first n 'linestocut'.
728
Usage: quickload(fname,linestocut=4)
729
Returns: array filled with data in fname
734
print fname,'read in.'
736
d = map(string.split,d)
737
print 'Done with string.split on lines.'
738
for i in range(len(d)):
739
d[i] = map(string.atoi,d[i])
740
print 'Conversion to ints done.'
743
def writedelimited(listoflists, delimiter, file, writetype='w'):
745
Writes a list of lists in columns, separated by character(s) delimiter
746
to specified file. File-overwrite is the default.
748
Usage: writedelimited(listoflists,delimiter,filename,writetype='w')
751
if type(listoflists[0]) not in [ListType,TupleType]:
752
listoflists = [listoflists]
753
outfile = open(file,writetype)
755
list2print = copy.deepcopy(listoflists)
756
for i in range(len(listoflists)):
757
if listoflists[i] == ['\n'] or listoflists[i]=='\n' or listoflists[i]=='dashes':
758
rowstokill = rowstokill + [i]
760
for row in rowstokill:
762
maxsize = [0]*len(list2print[0])
763
for row in listoflists:
764
if row == ['\n'] or row == '\n':
766
elif row == ['dashes'] or row == 'dashes':
767
dashes = [0]*len(maxsize)
768
for j in range(len(maxsize)):
770
outfile.write(pstat.linedelimited(dashes,delimiter))
772
outfile.write(pstat.linedelimited(row,delimiter))
777
def writecc(listoflists,file,writetype='w',extra=2):
779
Writes a list of lists to a file in columns, customized by the max
780
size of items within the columns (max size of items in col, +2 characters)
781
to specified file. File-overwrite is the default.
783
Usage: writecc(listoflists,file,writetype='w',extra=2)
786
if type(listoflists[0]) not in [ListType,TupleType]:
787
listoflists = [listoflists]
788
outfile = open(file,writetype)
790
list2print = copy.deepcopy(listoflists)
791
for i in range(len(listoflists)):
792
if listoflists[i] == ['\n'] or listoflists[i]=='\n' or listoflists[i]=='dashes':
793
rowstokill = rowstokill + [i]
795
for row in rowstokill:
797
maxsize = [0]*len(list2print[0])
798
for col in range(len(list2print[0])):
799
items = pstat.colex(list2print,col)
800
items = map(pstat.makestr,items)
801
maxsize[col] = max(map(len,items)) + extra
802
for row in listoflists:
803
if row == ['\n'] or row == '\n':
805
elif row == ['dashes'] or row == 'dashes':
806
dashes = [0]*len(maxsize)
807
for j in range(len(maxsize)):
808
dashes[j] = '-'*(maxsize[j]-2)
809
outfile.write(pstat.lineincustcols(dashes,maxsize))
811
outfile.write(pstat.lineincustcols(row,maxsize))
817
def writefc(listoflists,colsize,file,writetype='w'):
819
Writes a list of lists to a file in columns of fixed size. File-overwrite
822
Usage: writefc(listoflists,colsize,file,writetype='w')
825
if type(listoflists) == N.ndarray:
826
listoflists = listoflists.tolist()
827
if type(listoflists[0]) not in [ListType,TupleType]:
828
listoflists = [listoflists]
829
outfile = open(file,writetype)
831
list2print = copy.deepcopy(listoflists)
832
for i in range(len(listoflists)):
833
if listoflists[i] == ['\n'] or listoflists[i]=='\n' or listoflists[i]=='dashes':
834
rowstokill = rowstokill + [i]
836
for row in rowstokill:
838
n = [0]*len(list2print[0])
839
for row in listoflists:
840
if row == ['\n'] or row == '\n':
842
elif row == ['dashes'] or row == 'dashes':
844
for j in range(len(n)):
845
dashes[j] = '-'*(colsize)
846
outfile.write(pstat.lineincols(dashes,colsize))
848
outfile.write(pstat.lineincols(row,colsize))
854
def load(fname,lines_to_ignore=4,type='i'):
856
Load in huge, flat, 2D text files. Can handle differing line-lengths AND
857
can strip #/% on UNIX (or with a better NT grep). Requires wc, grep, and
858
mmapfile.lib/.pyd. Type can be 'i', 'f' or 'd', for ints, floats or doubles,
859
respectively. Lines_to_ignore determines how many lines at the start of the
860
file to ignore (required for non-working grep).
862
Usage: load(fname,lines_to_ignore=4,type='i')
863
Returns: numpy array of specified type
865
start = time.time() ## START TIMER
868
elif type in ['f','d']:
871
raise ValueError, "type can be 'i', 'f' or 'd' in load()"
873
## STRIP OUT % AND # LINES
874
tmpname = tempfile.mktemp()
875
if sys.platform == 'win32':
876
# NT VERSION OF GREP DOESN'T DO THE STRIPPING ... SIGH
877
cmd = "grep.exe -v \'%\' "+fname+" > "+tmpname
881
# UNIX SIDE SHOULD WORK
882
cmd = "cat "+fname+" | grep -v \'%\' |grep -v \'#\' > "+tmpname
886
## GET NUMBER OF ROWS, COLUMNS AND LINE-LENGTH, USING WC
887
wc = string.split(os.popen("wc "+tmpname).read())
888
numlines = int(wc[0]) - lines_to_ignore
890
if lines_to_ignore <> 0:
891
for i in range(lines_to_ignore):
892
junk = tfp.readline()
893
numcols = len(string.split(tfp.readline())) #int(float(wc[1])/numlines)
896
## PREPARE INPUT SPACE
897
a = N.zeros((numlines*numcols), type)
898
block = 65536 # chunk to read, in bytes
899
data = mmapfile.mmapfile(tmpname, '', 0)
900
if lines_to_ignore <> 0 and sys.platform == 'win32':
901
for i in range(lines_to_ignore):
902
junk = data.readline()
907
d = carryover + data.read(block)
908
cutindex = string.rfind(d,'\n')
909
carryover = d[cutindex+1:]
911
d = map(intype,string.split(d))
915
print "%d sec" % round(end-start,2)
918
return N.reshape(a,[numlines,numcols])
921
def find_dirs(sourcedir):
922
'''Finds and returns all directories in sourcedir
924
Usage: find_dirs(sourcedir)
925
Returns: list of directory names (potentially empty)
927
files = os.listdir(sourcedir)
930
if os.path.isdir(os.path.join(sourcedir,fname)):
940
def binget(fname,btype=None):
942
Loads a binary file from disk. Assumes associated hdr file is in same
943
location. You can force an unpacking type, or else it tries to figure
944
it out from the filename (4th-to-last character). Hence, readable file
947
1bin=int8, sbin=int16, ibin=int32, fbin=float32, dbin=float64, etc.
949
Usage: binget(fname,btype=None)
950
Returns: data in file fname of type btype
952
file = open(fname,'rb')
956
# if none given, assume character preceeding 'bin' is the unpacktype
960
bdata = N.fromstring(bdata,btype)
962
raise ValueError, "Bad unpacking type."
964
# force the data on disk to be LittleEndian (for more efficient PC/Linux use)
965
if not N.little_endian:
966
bdata = bdata.byteswap()
969
header = fname[:-3]+'hdr'
970
vals = get(header,0) # '0' means no missing-file warning msg
972
if type(vals[0]) == ListType: # it's an extended header
973
xsize = int(vals[0][0])
974
ysize = int(vals[0][1])
975
numslices = int(vals[0][2])
979
print "No (or bad) header file. Returning unshaped array."
980
return N.array(bdata)
984
def binput(outarray,fname,packtype=None,writetype='wb'):
986
Unravels outarray and writes the data to a file, always in LittleEndian
987
format, along with a header file containing the original data shape. Default
988
is overwrite the destination file. Tries to figure out packtype from
989
4th-to-last character in filename. Thus, the routine understands these
992
1bin=int8, sbin=int16, ibin=int32, fbin=float32, dbin=float64, etc.
994
Usage: binput(outarray,filename,packtype=None,writetype='wb')
999
# a speck of error checking
1000
if packtype == N.int16 and outarray.dtype.char == 'f':
1001
# check to see if there's data loss
1002
if max(N.ravel(outarray)) > 32767 or min(N.ravel(outarray))<-32768:
1003
print "*** WARNING: CONVERTING FLOAT DATA TO OUT-OF RANGE INT16 DATA"
1004
outdata = N.ravel(outarray).astype(packtype)
1006
# force the data on disk to be little_endian (for more efficient PC/Linux use)
1007
if not N.little_endian:
1008
outdata = outdata.byteswap()
1009
outdata = outdata.tostring()
1010
outfile = open(fname,writetype)
1011
outfile.write(outdata)
1014
# Now, write the header file
1016
suffixindex = string.rfind(fname,'.')
1017
hdrname = fname[0:suffixindex+2]+'hdr' # include .s or .f or .1 or whatever
1020
hdr = outarray.shape
1022
outfile = open(hdrname,'w')
1023
outfile.write(pstat.list2string(hdr))
1027
def getafniparam(headfilename,paramname):
1029
Loads in an AFNI header file, and returns the values of 'paramname'.
1031
Usage: getafniparam(headfile,paramname)
1032
Returns: appropriate "type" for params, or None if fails
1034
if headfilename[-4:] == 'BRIK': # if asked for BRIK, change it to HEAD
1035
headfilename = headfilename[:-4]+'HEAD'
1036
d = get(headfilename)
1037
lines = open(headfilename,'r').readlines()
1038
for i in range(len(lines)):
1039
if string.find(lines[i],paramname) <> -1:
1043
for j in range(i+2,len(lines)):
1044
for k in range(len(d[j])):
1045
if type(d[j][k]) == StringType:
1046
result = d[j][k][1:count]
1049
result.append(d[j][k])
1057
def add2afnihistory(headfilename,newtext):
1059
Adds 'newtext' to HISTORY_NOTE in afni file specified in headfilename.
1061
Usage: add2afnihistory(headfile,newtext)
1064
if headfilename[-4:] == 'BRIK': # if asked for BRIK, change it to HEAD
1065
headfilename = headfilename[:-4]+'HEAD'
1066
d = get(headfilename)
1067
lines = open(headfilename,'r').readlines()
1068
for i in range(len(lines)):
1069
if string.find(lines[i],'HISTORY_NOTE') <> -1:
1070
bytecount = d[i+1][-1]
1071
oldstr = lines[i+2][:-2]
1072
date = '[python:*** %s] ' %time.asctime()
1073
lines[i+2] = oldstr +'\\n' +date +newtext +'~\n'
1074
lines[i+1] = ' count = %s\n' %str(len(lines[i+2]))
1075
f = open(headfilename,'w')
1081
def array2afni(d,brikprefix,voltype=None,TR=2000,sliceorder='seqplus',geomparent=None,view=None,corrlength=1,briklabels=None,historytext=None):
1083
Converts an array 'd' to an AFNI BRIK/HEAD combo via putbin and to3d. Tries to
1084
guess the AFNI volume type
1086
voltype = {'-anat','-epan','-fim'}
1087
geomparent = filename of the afni BRIK file with the same geometry
1088
view = {'tlrc', 'acpc' or 'orig'}
1089
corrlength = # of images used in the (single-timeseries) correlation (for fico)
1090
briklabels = list of names (strings) to use for brick labels
1091
historytext = string to be appended to the history file, if any
1093
Usage: array2afni(d,brikprefix,voltype=None,TR=2000,
1094
sliceorder='seqplus',geomparent=None,view=None,
1095
corrlength=1,briklabels=None,historytext=None)
1098
# converts numpy typecode()s into appropriate strings for to3d command line
1099
typecodemapping = {'c':'b', # character
1100
'B':'b', # UnsignedInt8
1101
'f':'f', # float0, float8, float16, float32
1103
'b':'b', # int0, int8
1108
# Verify that the data is proper size (3- or 4-D)
1109
if len(d.shape) not in [3,4]:
1110
raise ValueError, "A 3D or 4D array is required for array2afni() ... %s" %d.shape
1112
# Save out the array to a binary file, homebrew style
1113
if d.dtype.char == N.float64:
1116
outcode = d.dtype.char
1117
tmpoutname = 'afnitmp.%sbin' % outcode
1118
binput(d.astype(outcode),tmpoutname)
1120
if len(d.shape) == 3: # either anatomy or functional
1121
if d.dtype.char in ['s','i','l']: # if floats, assume functional
1125
else: # 4D dataset, must be anatomical timeseries (epan)
1127
if voltype[0] != '-':
1128
voltype = '-'+voltype
1129
if len(d.shape) == 3: # either anatomy or functional
1133
elif len(d.shape) == 4:
1134
if voltype=='-fico':
1136
d = N.reshape(d,[d.shape[0]*d.shape[1],d.shape[2],d.shape[3]])
1138
timestr = '-statpar %s 1 1 ' % corrlength
1140
timepts = d.shape[0]
1142
timestr = '-time:zt %d %d %0.3f %s ' % (slices,timepts,TR,sliceorder)
1144
cmd = 'to3d %s -prefix %s -session . ' % (voltype, brikprefix)
1147
cmd += '-view %s ' % view
1149
cmd += '-geomparent %s ' % geomparent
1151
cmd += '3D%s:0:0:%d:%d:%d:%s' % (typecodemapping[d.dtype.char],d.shape[-1],d.shape[-2],slices*timepts,tmpoutname)
1154
os.remove(tmpoutname)
1155
os.remove(tmpoutname[:-3]+'hdr')
1157
if len(d.shape)==4 and briklabels:
1159
for label in briklabels:
1160
names += str(label)+'~'
1162
appendstr = '''\n\ntype = string-attribute
1165
'%s''' % (count, names)
1166
f = open('%s+%s.HEAD' %(brikprefix,view), 'a')
1171
add2afnihistory('%s+%s.HEAD'%(brikprefix,view),historytext)