3
// hack to create voxbo glm's
4
// Copyright (c) 1998-2010 by The VoxBo Development Team
6
// VoxBo is free software: you can redistribute it and/or modify it
7
// under the terms of the GNU General Public License as published by
8
// the Free Software Foundation, either version 3 of the License, or
9
// (at your option) any later version.
11
// VoxBo is distributed in the hope that it will be useful, but
12
// WITHOUT ANY WARRANTY; without even the implied warranty of
13
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14
// General Public License for more details.
16
// You should have received a copy of the GNU General Public License
17
// along with VoxBo. If not, see <http://www.gnu.org/licenses/>.
19
// For general information on VoxBo, including the latest complete
20
// source code and binary distributions, manual, and associated files,
21
// see the VoxBo home page at: http://www.voxbo.org/
23
// original version written by Dan Kimberg
24
// written to accomplish the same thing as some original IDL code by
35
#include "vbjobspec.h"
44
const int BUFLEN=1024;
45
bool f_ignorewarnings,f_validateonly,f_setuponly;
51
string name; // name of the sequence
52
string dirname; // root directory
53
string stem; // stem name for analysis
54
vector<string>scanlist; // list of scans (data files)
55
int lows,highs; // low and high frequencies to remove
56
string middles; // middle frequencies to remove
57
int pieces; // n pieces for matrix jobs
58
string kernelname; // temporal smoothing kernel
59
double kerneltr; // sampling interval of kernel file
60
string noisemodel; // noise model
61
string refname; // file containing reference function
62
string glmfile; // filename of .glm file
63
string gmatrix; // g matrix file
64
string email; // email address
65
vector<string> contrasts; // list of contrasts to add to contrasts.txt
67
bool auditflag,meannorm,emailflag,meannormset,driftcorrect,driftcorrectset;
71
bool rfxgflag; // should we just dummy up a G matrix with all 1's?
76
void FixRelativePaths();
78
void CreateGLMJobs2();
79
vector<string> CreateGLMScript();
80
int WriteGLMFile(string filename="");
81
int WriteAndSubmitJobs();
83
int createsamplefiles();
88
GLMParams global,local,*gp;
90
GLMConfig(string fname,double tr,double kr);
91
int ParseFile(string filename);
94
void vbmakeglm_help();
95
void vbmakeglm_sample(FILE *);
96
void vbmakeglm_make_sample(const string &fname);
97
void vbmakeglm_show_sample();
100
main(int argc,char **argv)
107
args.Transfer(argc-1,argv+1);
108
if (args.size() == 0) {
113
vbmakeglm_show_sample();
118
vbmakeglm_make_sample(args[1]);
120
printf("vbmakeglm -x needs a filename to store the sample glm file\n");
129
vector<string> filelist;
131
double x_tr=0.0,x_kr=0.0;
133
for (i=0; i<args.size(); i++) {
138
else if (args[i]=="-v") {
142
else if (args[i]=="-s") {
146
else if (args[i]=="--run") {
148
if (f_run<1) f_run=1;
150
else if (args[i].compare(0,6,"--run=")==0) {
151
f_run=strtol(args[i].substr(6));
152
if (f_run<1) f_run=1;
154
else if (args[i]=="-n" && (i+1)<args.size()) {
159
else if (args[i]=="-tr" && (i+1)<args.size()) {
160
x_tr=strtod(args[i+1]);
164
else if (args[i]=="-kr" && (i+1)<args.size()) {
165
x_kr=strtod(args[i+1]);
170
filelist.push_back(args[i]);
173
for (i=0; i<(int)filelist.size(); i++) {
174
GLMConfig *glmc=new GLMConfig(filelist[i],x_tr,x_kr);
180
// this is the sole constructor, it chews on the passed arguments
182
GLMConfig::GLMConfig(string fname,double tr,double kr)
187
global.glmfile=fname;
221
GLMParams::Validate()
233
printf("[I] vbmakeglm: validating GLM %s...\n",name.c_str());
235
if (!(dirname.size())) {
236
printf("[E] vbmakeglm: no dirname specified for glm %s\n",name.c_str());
240
for (i=0; i<(int)scanlist.size(); i++) {
242
mytes->ReadHeader(scanlist[i].c_str());
243
if (!(mytes->header_valid)) {
244
printf("[E] vbmakeglm: couldn't get info on data file %s.\n",scanlist[i].c_str());
247
else if ((string)"tes1" != mytes->GetFileFormat().getSignature()) {
248
printf("[E] vbmakeglm: couldn't read 4D data from file %s.\n",scanlist[i].c_str());
251
else if (mytes->GetFileFormat().getDimensions() != 4) {
252
printf("[E] vbmakeglm: %s isn't a 4D file.\n",scanlist[i].c_str());
255
if (mytes->header_valid) {
256
timecount += mytes->dimt;
257
if (i==0 && TR<1.0) {
258
if (mytes->voxsize[3]>0.0)
259
TR=mytes->voxsize[3];
264
// the number of time points ought to be the same as orderg
268
if (lows < 0 || lows > (orderg/2)) {
269
printf("[W] vbmakeglm: removing %d low frequencies is a little suspect\n",lows);
273
if (highs < 0 || highs > 20) {
274
printf("[W] vbmakeglm: removing %d high frequencies is a little suspect\n",highs);
278
if (noisemodel.size()) {
279
Vec p(noisemodel.c_str());
281
printf("[E] vbmakeglm: your 1/f parameter file %s has the wrong number of elements\n",xfilename(noisemodel).c_str());
286
if (kernelname.size()) {
288
if (IRF.ReadFile(kernelname.c_str())) {
290
tmps << "vbmakeglm: your HRF kernel file " << kernelname << " doesn't exist.";
291
printErrorMsg(VB_ERROR,tmps.str());
294
else if (IRF.size() < 6 || IRF.size() > 20) {
296
tmps << "vbmakeglm: your HRF kernel file " << kernelname << " has "
297
<< IRF.size() << " elements, which seems a little suspect.";
298
printErrorMsg(VB_WARNING,tmps.str());
301
else if (kerneltr < 1.0) {
302
string trline=GetHeader(IRF.header,"TR(msecs)");
305
tmpl.ParseLine(trline);
307
kerneltr=strtod(tmpl[1]);
312
// load the G matrix, check its order and see if it has an intercept of interest
313
if (gmatrix.size()) {
314
VBMatrix gmat(gmatrix);
315
// gmat.MakeInCore();
316
if (gmat.m <=0 || gmat.n <= 0) {
318
tmps << "vbmakeglm: couldn't read G (design) matrix " << gmatrix << ".";
319
printErrorMsg(VB_ERROR,tmps.str());
324
if (gmat.n==1 && meannorm) {
326
tmps << "vbmakeglm: you have a single covariate and the mean norm flag set - make sure you don't "
327
<< "mean normalize data for second tier (group rfx) analyses.";
328
printErrorMsg(VB_WARNING,tmps.str());
333
if (rfxgflag && meannorm) {
335
tmps << "vbmakeglm: you have the makerrandfxg flag and the mean norm flag set - make sure you don't "
336
<< "mean normalize data for second tier (group rfx) analyses.";
337
printErrorMsg(VB_WARNING,tmps.str());
341
if (meannorm && noisemodel.size()==0) {
343
tmps << "vbmakeglm: you have the mean norm flag set, and no noise model - make sure you really want"
344
<< "to mean normalize these data.";
345
printErrorMsg(VB_WARNING,tmps.str());
351
tmps << "vbmakeglm: no meannorm flag set -- defaulting to no mean normalization";
352
printErrorMsg(VB_WARNING,tmps.str());
355
if (orderg != timecount) {
356
printf("[E] vbmakeglm: orderg (%d) doesn't match the number of timepoints (%d)\n",
361
// FIXME - should actually check order of G matrix
362
if (gmatrix.size()==0 && !rfxgflag) {
363
printf("[E] vbmakeglm: no valid G (design) matrix specified\n");
369
if (warns && !f_ignorewarnings)
374
tmps << "vbmakeglm: GLM " << name << " is good to go.";
375
printErrorMsg(VB_INFO,tmps.str());
380
tmps << "vbmakeglm: TR not set, using default of 2000ms.";
381
printErrorMsg(VB_WARNING,tmps.str());
383
if (kerneltr<1.0 && kernelname.size()) {
386
tmps << "vbmakeglm: HRF TR (sampling rate) not set, using default of 2000ms.";
387
printErrorMsg(VB_WARNING,tmps.str());
389
// the following added to double-check on the TR bugs
390
printf("[I] vbmakeglm: your data TR is %f\n",TR);
391
printf("[I] vbmakeglm: your HRF kernel TR (sampling rate) is %f\n",kerneltr);
395
GLMParams::CreateGLMDir()
397
string fname; // temporary variable for wherever we need to build a filename
403
stemstuff.SetSeparator(" \t\n\r/");
404
stemstuff.ParseLine(dirname.c_str());
405
if (!stemstuff.size()) {
406
tmps << "vbmakeglm: couldn't find parent of GLM directory" << dirname << ".";
407
printErrorMsg(VB_ERROR,tmps.str());
410
stem=dirname+(string)"/"+stemstuff[stemstuff.size()-1];
413
mkdir(dirname.c_str(),0777);
414
if (stat(dirname.c_str(),&st)) {
415
tmps << "vbmakeglm: couldn't get to GLM directory " << dirname << ".";
416
printErrorMsg(VB_ERROR,tmps.str());
421
ofstream subfile(fname.c_str());
423
tmps << "vbmakeglm: couldn't create sub file " << fname << ".";
424
printErrorMsg(VB_ERROR,tmps.str());
427
subfile << ";VB98\n;TXT1\n;\n; Subject list generated by vbmakeglm\n;\n\n";
428
for (i=0; i<(int)scanlist.size(); i++)
429
subfile << scanlist[i] << endl;
433
if (refname.size()) {
434
if (copyfile(refname,stem+".ref")) {
435
tmps << "vbmakeglm: couldn't copy reference function " << refname << ".";
436
printErrorMsg(VB_ERROR,tmps.str());
441
// create glm file for potential future use
442
if (glmfile.size()) {
443
if (WriteGLMFile(stem+".glm")) {
444
tmps << "vbmakeglm: couldn't write glm configuration file " << (stem+".glm") << ".";
445
printErrorMsg(VB_WARNING,tmps.str());
450
// copy G if necessary
451
if (gmatrix.size() > 0) {
452
if (copyfile(gmatrix,stem+".G")) {
453
printf("[E] vbmakeglm: couldn't copy G matrix file %s.\n",gmatrix.c_str());
459
ofstream gstr(gmatrix.c_str(),ios::binary);
461
gstr << "VB98\nMAT1\n";
462
gstr << "DataType:\tFloat\n";
463
gstr << "VoxDims(XY):\t1\t" << orderg << endl << endl;
464
gstr << "# This G matrix generated automatically by vbmakeglm\n\n";
465
gstr << "Parameter:\t0\tInterest\tEffect\n";
468
for (i=0; i<orderg; i++)
470
if (my_endian() != ENDIAN_BIG)
472
for (i=0; i<(int)(orderg * sizeof(float)); i++) {
473
gstr << *((unsigned char *)pts+i);
478
createsamplefiles(); // after G matrix, so we can count vars of interest
479
return 0; // no error!
483
GLMParams::createsamplefiles()
485
// first find variables of interest
488
glmi.getcovariatenames();
489
string fname=dirname+"/contrasts.txt";
490
vector<string> interestnames;
491
if (access(fname.c_str(),R_OK) || contrasts.size()) {
492
ofstream ofile(fname.c_str());
494
ofile << "# contrasts.txt\n";
495
ofile << "# in this file you can define contrasts among your covariates of interest\n";
496
if (glmi.cnames.size()) {
497
ofile << "# your covariates of interest are:\n";
498
for (size_t i=0; i<glmi.cnames.size(); i++) {
499
if (glmi.cnames[i][0]=='I') {
500
ofile << "# "<<strnum(i)<<": "<<glmi.cnames[i].c_str()+1<<endl;
501
interestnames.push_back(glmi.cnames[i].substr(1));
505
ofile << "# you can specify a complete contrast as follows:\n#\n";
506
ofile << "# <name> <scale> vec"; // 0 0 1 0\n";
508
for (size_t i=1; i<interestnames.size(); i++)
510
ofile << endl << "#\n";
511
ofile << "# (with one value for each covariate of interest)\n";
513
ofile << "# lines beginning with a '#' are comments!\n";
515
ofile << "# the following simple contrasts are provided for your convenience:\n";
517
for (size_t i=0; i<interestnames.size(); i++) {
518
ofile << interestnames[i] << " t vec";
519
for (size_t j=0; j<interestnames.size(); j++) {
520
if (j==i) ofile << " 1";
525
if (contrasts.size()) {
526
ofile << "\n# the following contrasts were specified:\n";
528
for (size_t i=0; i<contrasts.size(); i++) {
529
if (glmi.parsecontrast(contrasts[i]))
530
printf("[W] vbgmakeglm: couldn't parse contrast: %s\n",contrasts[i].c_str());
532
ofile << contrasts[i] << endl;
539
fname=dirname+"/averages.txt";
540
if (access(fname.c_str(),R_OK)) {
541
ofstream ofile(fname.c_str());
543
ofile << "# averages.txt\n";
545
ofile << "# In this file you can specify one or more ways to trial-average your data.\n";
546
ofile << "# You can also block-average or whatever else you need, we just call it\n";
547
ofile << "# trial averaging generically.\n";
549
ofile << "# Each trial average needs a separate section that looks like the following:\n";
551
ofile << "# average <name>\n";
552
ofile << "# units <time/vols>\n";
553
ofile << "# interval <ms/vols>\n";
554
ofile << "# nsamples <n>\n";
555
ofile << "# tr <ms>\n";
556
ofile << "# trial <n>...\n";
557
ofile << "# trialset <first> <interval> <count>\n";
560
ofile << "# Here's what these parameters mean:\n";
562
ofile << "# units: whether the other parameters are in volumes or seconds\n";
563
ofile << "# interval: interval in time or volumes between samples within the trial\n";
564
ofile << "# nsamples: number of time points to include per trial\n";
565
ofile << "# tr: sets the TR if you're using time units\n";
567
ofile << "# The remaining options are two ways to indicate when the trials begin.\n";
568
ofile << "# If your trials are evenly spaced, use 'trialset,' otherwise use 'trial'\n";
570
ofile << "# trialset: specify the start of the first trial, the interval between trial\n";
571
ofile << "# onsets, and the trial count\n";
572
ofile << "# trial: each trial line lists one or more start times/vols for a trial\n";
573
ofile << "# (you can include multiple trial lines to help you keep the file neat)\n";
575
ofile << "# --> for trial and trialsets, the first volume is volume 0 (also time 0)\n";
576
ofile << "# --> both time and volume values can be floating point\n";
578
ofile << "# Total data points for this GLM: "<<orderg<<endl;
579
ofile << "# Your TR in ms: "<<TR<<endl;
589
GLMConfig::ParseFile(string filename)
596
infile.open(filename.c_str());
598
tmps << "vbmakeglm: couldn't open " << filename;
599
printErrorMsg(VB_ERROR,tmps.str());
602
while (infile.getline(buf,BUFLEN,'\n')) {
605
if (args[0][0]=='#' || args[0][0] == ';')
607
if (args[0]=="name" || args[0]=="glm") {
608
local=global; // start with a copy of the global
610
gp->name=args.Tail();
612
else if (args[0]=="dirname") {
613
if (args.size() == 2)
617
tmps << "vbmakeglm: dirname takes exactly one argument";
618
printErrorMsg(VB_ERROR,tmps.str());
621
else if (args[0]=="scan") {
622
if (args.size() != 2) {
624
tmps << "vbmakeglm: scan takes exactly one argument";
625
printErrorMsg(VB_ERROR,tmps.str());
628
gp->scanlist.push_back(args[1]);
630
else if (args[0]=="lows") {
631
if (args.size() == 2)
632
gp->lows=strtol(args[1]);
635
tmps << "vbmakeglm: lows takes exactly one argument";
636
printErrorMsg(VB_ERROR,tmps.str());
639
else if (args[0]=="highs") {
640
if (args.size() == 2)
641
gp->highs=strtol(args[1]);
644
tmps << "vbmakeglm: highs takes exactly one argument";
645
printErrorMsg(VB_ERROR,tmps.str());
648
else if (args[0]=="middles") {
652
else if (args[0]=="pieces") {
653
if (args.size() == 2)
654
gp->pieces=strtol(args[1]);
656
printf("[E] vbmakeglm: pieces takes exactly one argument\n");
659
else if (args[0]=="orderg") {
660
if (args.size() == 2)
661
gp->orderg=strtol(args[1]);
664
tmps << "vbmakeglm: orderg takes exactly one argument";
665
printErrorMsg(VB_ERROR,tmps.str());
668
else if (args[0]=="refname") {
669
if (args.size() == 2)
672
printf("[E] vbmakeglm: reference takes exactly one argument\n");
673
printErrorMsg(VB_ERROR,tmps.str());
676
else if (args[0]=="kernel") {
677
if (args.size() == 2)
678
gp->kernelname=args[1];
679
else if (args.size()==3) {
680
gp->kernelname=args[1];
681
gp->kerneltr=strtod(args[2]);
685
tmps << "vbmakeglm: kernel takes one or two arguments";
686
printErrorMsg(VB_ERROR,tmps.str());
689
else if (args[0]=="noisemodel") {
690
if (args.size() == 2)
691
gp->noisemodel=args[1];
694
tmps << "vbmakeglm: noisemodel takes exactly one argument";
695
printErrorMsg(VB_ERROR,tmps.str());
698
else if (args[0]=="gmatrix") {
699
if (args.size() == 2)
703
tmps << "vbmakeglm: gmatrix takes exactly one argument";
704
printErrorMsg(VB_ERROR,tmps.str());
707
else if (args[0]=="makerandfxg") {
709
if (args.size() > 1) {
711
tmps << "vbmakeglm: arguments to makerandfxg ignored";
712
printErrorMsg(VB_WARNING,tmps.str());
715
else if (args[0]=="tr") {
716
if (args.size() == 2)
717
gp->TR=strtol(args[1]);
720
tmps << "vbmakeglm: tr takes exactly one argument";
721
printErrorMsg(VB_ERROR,tmps.str());
724
else if (args[0]=="pri") {
725
if (args.size() == 2)
726
gp->pri=strtol(args[1]);
729
tmps << "vbmakeglm: pri takes exactly one argument";
730
printErrorMsg(VB_ERROR,tmps.str());
733
else if (args[0]=="audit") {
734
if (args.size() != 2) {
736
tmps << "vbmakeglm: audit takes exactly one argument";
737
printErrorMsg(VB_ERROR,tmps.str());
739
if (args[1]=="yes" || args[1]=="1")
741
else if (args[1]=="no" || args[1]=="0")
745
tmps << "vbmakeglm: unrecognized value for audit flag (should be yes/no)";
746
printErrorMsg(VB_ERROR,tmps.str());
749
else if (args[0]=="meannorm") {
750
if (args.size() != 2) {
752
tmps << "vbmakeglm: meannorm takes exactly one argument";
753
printErrorMsg(VB_ERROR,tmps.str());
755
if (args[1]=="yes" || args[1]=="1") {
759
else if (args[1]=="no" || args[1]=="0") {
764
else if (args[0]=="driftcorrect") {
765
if (args.size() != 2) {
767
tmps << "vbmakeglm: driftcorrect takes exactly one argument";
768
printErrorMsg(VB_ERROR,tmps.str());
770
if (args[1]=="yes" || args[1]=="1") {
772
gp->driftcorrectset=1;
774
else if (args[1]=="no" || args[1]=="0") {
776
gp->driftcorrectset=1;
780
tmps << "vbmakeglm: unrecognized value for driftcorrect flag (should be yes/no)";
781
printErrorMsg(VB_ERROR,tmps.str());
784
else if (args[0]=="email") {
786
if (args.size() == 2)
789
else if (args[0]=="go" || args[0]=="end") {
790
if (!glmname.size() || glmname==gp->name) {
791
gp->FixRelativePaths();
793
if (gp->valid && !f_validateonly) {
794
if (!gp->CreateGLMDir()) {
798
gp->CreateGLMJobs2();
799
runseq(vbp,gp->seq,f_run);
803
gp->WriteAndSubmitJobs();
810
else if (args[0]=="contrast") {
811
gp->contrasts.push_back(args.Tail());
813
else if (args[0]=="endx" || args[0]=="noend") {
816
else if (args[0]=="include") {
817
if (args.size() == 2)
822
tmps << "vbmakeglm: unrecognized keyword " << args[0];
823
printErrorMsg(VB_ERROR,tmps.str());
828
return 0; // no error!
832
GLMParams::CreateGLMJobs()
835
int rowstart,rowfinish,mergeflag;
841
// set pieces heuristically if the user didn't
843
// FIXME is this any good?
844
int cellcount=600000;
845
pieces=(int)ceil((double)orderg*orderg/cellcount);
846
if (pieces > orderg) pieces=orderg;
847
if (pieces<1) pieces=1;
849
int rowsperjob=orderg/pieces;
850
// are we going to need the merge jobs?
851
if (rowsperjob < orderg)
858
js.jobtype="vb_makefilter";
859
js.arguments["outfile"]=stem+".ExoFilt";
860
// sprintf(tmp,"lowflag -lf %d",lows);
861
js.arguments["lowflag"]=(string)"-lf "+strnum(lows);
863
js.arguments["middleflag"]=(string)"-mf "+middles;
865
js.arguments["middleflag"]="";
866
sprintf(tmp,"highflag -hf %d",highs);
867
js.arguments["highflag"]=(string)"-hf "+strnum(highs);
868
if (kernelname.size())
869
js.arguments["kernelflag"]=(string)"-k "+kernelname;
871
js.arguments["kernelflag"]="";
872
js.arguments["ndata"]=strnum(orderg);
873
sprintf(tmp,"tr %f",TR);
874
js.arguments["tr"]=strnum(TR);
875
js.magnitude=0; // set it to something!
876
js.name="make exofilt";
878
int n_exofilt=js.jnum;
881
// make the noisemodel
883
js.jobtype="vb_makenoisemodel";
884
js.arguments["outfile"]=stem+".IntrinCor";
885
if (noisemodel.size())
886
js.arguments["noisemodelflag"]="-n "+noisemodel;
887
js.arguments["ndata"]=strnum(orderg);
888
js.arguments["tr"]=strnum(TR);
889
js.magnitude=0; // set it to something!
890
js.name="make noisemodel";
892
int n_noisemodel=js.jnum;
895
// make the KG matrix
897
js.jobtype="makematkg";
898
js.arguments["stem"]=stem;
899
js.magnitude=0; // set it to something!
901
js.waitfor.insert(n_exofilt);
906
// make the F1 matrix
908
js.jobtype="matpinv";
909
js.arguments["in"]=stem+".KG";
910
js.arguments["out"]=stem+".F1";
911
js.magnitude=0; // set it to something!
913
js.waitfor.insert(n_kg);
918
// make the R matrix = I-(KG)(F1)
920
js.jobtype="matimxy";
921
js.arguments["in1"]=stem+".KG";
922
js.arguments["in2"]=stem+".F1";
923
js.arguments["out"]=stem+".R";
924
js.magnitude=0; // set it to something!
926
js.waitfor.insert(n_f1);
933
js.jobtype="makematk";
934
js.arguments["stem"]=stem;
935
js.magnitude=0; // set it to something!
937
js.waitfor.insert(n_exofilt);
938
js.waitfor.insert(n_noisemodel);
943
// create an empty V matrix
945
js.jobtype="matzeros";
946
js.arguments["name"]=stem+".V";
947
js.arguments["cols"]=strnum(orderg);
948
js.arguments["rows"]=strnum(orderg);
949
js.magnitude=0; // set it to something!
950
js.name="GLM-Vcreate";
952
int n_vcreate=js.jnum;
955
// build V matrix (V=KKt)
956
set<int32> n_vpieces;
961
rowfinish=rowstart+rowsperjob;
962
if (rowfinish > orderg-1)
966
js.arguments["in"]=stem+".K";
967
js.arguments["out"]=stem+".V";
968
js.arguments["col1"]=strnum(rowstart);
969
js.arguments["col2"]=strnum(rowfinish);
970
js.magnitude=0; // set it to something!
971
sprintf(tmp,"GLM-V%d",i++);
973
js.waitfor.insert(n_k);
974
js.waitfor.insert(n_vcreate);
976
n_vpieces.insert(js.jnum);
977
n_v=js.jnum; // in case it's just one piece
980
rowstart+=rowsperjob+1;
981
} while (rowstart < orderg);
986
js.jobtype="matmerge";
987
js.arguments["out"]=stem+".V";
988
js.name="GLM-MergeV";
989
js.waitfor=n_vpieces;
998
js.arguments["vmatrix"]=stem+".V";
999
js.arguments["kgmatrix"]=stem+".KG";
1000
js.arguments["outfile"]=stem+".F3";
1001
js.magnitude=0; // set it to something!
1002
js.name="GLM-makeF3";
1003
js.waitfor.insert(n_kg);
1004
js.waitfor.insert(n_v);
1009
// create an empty RV matrix
1011
js.jobtype="matzeros";
1012
js.arguments["name"]=stem+".RV";
1013
js.arguments["cols"]=strnum(orderg);
1014
js.arguments["rows"]=strnum(orderg);
1015
js.magnitude=0; // set it to something!
1016
js.name="GLM-RVcreate";
1018
int n_rvcreate=js.jnum;
1022
set<int32> n_rvpieces;
1026
rowfinish=rowstart+rowsperjob;
1027
if (rowfinish > orderg-1)
1031
js.arguments["in1"]=stem+".R";
1032
js.arguments["in2"]=stem+".V";
1033
js.arguments["out"]=stem+".RV";
1034
js.arguments["col1"]=strnum(rowstart);
1035
js.arguments["col2"]=strnum(rowfinish);
1036
js.magnitude=0; // set it to something!
1037
sprintf(tmp,"GLM-RV%d",i);
1039
js.waitfor.insert(n_r);
1040
js.waitfor.insert(n_v);
1041
js.waitfor.insert(n_rvcreate);
1045
n_rvpieces.insert(js.jnum);
1046
n_rv=js.jnum; // in case it's one piece
1048
rowstart+=rowsperjob+1;
1049
} while (rowstart < orderg);
1054
js.jobtype="matmerge";
1055
js.arguments["out"]=stem+".RV";
1056
js.name="GLM-MergeRV";
1057
js.waitfor=n_rvpieces;
1063
// create an empty RVRV matrix
1065
js.jobtype="matzeros";
1066
js.arguments["name"]=stem+".RVRV";
1067
js.arguments["cols"]=strnum(orderg);
1068
js.arguments["rows"]=strnum(orderg);
1069
js.magnitude=0; // set it to something!
1070
js.name="GLM-RVRVcreate";
1072
int n_rvrvcreate=js.jnum;
1076
set<int32> n_rvrvpieces;
1080
rowfinish=rowstart+rowsperjob;
1081
if (rowfinish > orderg-1)
1085
js.arguments["in"]=stem+".RV";
1086
js.arguments["out"]=stem+".RVRV";
1087
js.arguments["col1"]=strnum(rowstart);
1088
js.arguments["col2"]=strnum(rowfinish);
1089
js.magnitude=0; // set it to something!
1090
sprintf(tmp,"GLM-RVRV%d",i);
1092
js.waitfor.insert(n_rvrvcreate);
1093
js.waitfor.insert(n_rv);
1095
n_rvrvpieces.insert(js.jnum);
1099
rowstart+=rowsperjob+1;
1100
} while (rowstart < orderg);
1105
js.jobtype="matmerge";
1106
js.arguments["out"]=stem+".RVRV";
1107
js.name="GLM-MergeRVRV";
1108
js.waitfor=n_rvrvpieces;
1116
js.jobtype="comptraces";
1117
js.arguments["stem"]=stem;
1118
js.magnitude=0; // set it to something!
1119
js.name="GLM-traces";
1120
js.waitfor.insert(n_rvrv);
1122
int n_traces=js.jnum;
1125
// tes regression steps
1126
set<int32> n_regstep;
1127
for (i=0; i<pieces; i++) {
1129
js.jobtype="vbregress";
1130
js.arguments["stem"]=stem;
1131
js.arguments["steps"]=strnum(pieces);
1132
js.arguments["index"]=strnum(i+1);
1136
if (this->driftcorrect)
1138
js.arguments["flags"]=flags;
1139
js.magnitude=0; // set it to something!
1140
sprintf(tmp,"regress(%d/%d)",i+1,pieces);
1142
js.waitfor.insert(n_traces);
1143
js.waitfor.insert(n_f1);
1144
js.waitfor.insert(n_r);
1145
js.waitfor.insert(n_exofilt);
1147
n_regstep.insert(js.jnum);
1154
js.jobtype="vbmerge4d";
1155
js.arguments["stem"]=stem;
1156
js.arguments["ext"]="prm";
1157
js.magnitude=0; // set it to something!
1158
js.name="mergeparams";
1159
js.waitfor=n_regstep;
1165
js.jobtype="vbmerge4d";
1166
js.arguments["stem"]=stem;
1167
js.arguments["ext"]="res";
1168
js.magnitude=0; // set it to something!
1169
js.name="mergeparams";
1170
js.waitfor=n_regstep;
1175
// make these two jobs stand in for regstep in the waitfor chain
1177
n_regstep.insert(jobnum-2);
1178
n_regstep.insert(jobnum-1);
1181
// create smoothness estimate
1184
js.arguments["stem"]=stem;
1185
js.magnitude=0; // set it to something!
1187
js.waitfor=n_regstep;
1197
js.jobtype="vb_auditglm";
1198
js.arguments["glmdir"]=xdirname(stem);
1199
js.magnitude=0; // set it to something!
1200
js.name="GLM-audit";
1201
js.waitfor=n_regstep;
1202
js.waitfor.insert(n_se);
1209
if (emailflag && email.size()) {
1211
js.jobtype="notify";
1212
js.arguments["email"]=email;
1213
js.arguments["msg"]="Your GLM has been solved.";
1214
js.magnitude=0; // set it to something!
1216
js.waitfor.insert(n_se);
1218
js.waitfor.insert(n_audit);
1219
js.waitfor.insert(n_f3);
1226
GLMParams::CreateGLMJobs2()
1229
int rowstart,rowfinish,mergeflag;
1230
char tmp[STRINGLEN];
1235
// set pieces heuristically if the user didn't
1237
// FIXME is this any good?
1238
int cellcount=600000;
1239
pieces=(int)ceil((double)orderg*orderg/cellcount);
1240
if (pieces > orderg) pieces=orderg;
1241
if (pieces<1) pieces=1;
1243
int rowsperjob=orderg/pieces;
1244
// are we going to need the merge jobs?
1245
if (rowsperjob < orderg)
1252
js.jobtype="shellcommand";
1253
js.arguments["command"]=str(format("vbmakefilter -e %s.ExoFilt -lf %d -hf %d %s %s -t %d %f")
1254
%stem%lows%highs%(middles.size()?"-mf "+middles:"")
1255
%(kernelname.size()?"-k "+kernelname:"")
1257
js.name="make exofilt";
1259
int n_exofilt=js.jnum;
1262
// make the noisemodel
1264
js.jobtype="shellcommand";
1265
js.arguments["command"]=str(format("vbmakefilter -i %s.IntrinCor %s -t %d %f")
1266
%stem%(noisemodel.size()?"-n "+noisemodel:"")%orderg%TR);
1267
js.name="make noisemodel";
1269
int n_noisemodel=js.jnum;
1272
// make the KG matrix
1274
js.jobtype="shellcommand";
1275
js.arguments["command"]=str(format("makematkg -m %s")%stem);
1277
js.waitfor.insert(n_exofilt);
1282
// make the F1 matrix
1284
js.jobtype="shellcommand";
1285
js.arguments["command"]=str(format("vbmm2 -pinv %1%.KG %1%.F1")%stem);
1287
js.waitfor.insert(n_kg);
1292
// make the R matrix = I-(KG)(F1)
1294
js.jobtype="shellcommand";
1295
js.arguments["command"]=str(format("vbmm2 -imxy %1%.KG %1%.F1 %1%.R")%stem);
1297
js.waitfor.insert(n_f1);
1302
// create a K matrix
1304
js.jobtype="shellcommand";
1305
js.arguments["command"]=str(format("makematk -m %s")%stem);
1307
js.waitfor.insert(n_exofilt);
1308
js.waitfor.insert(n_noisemodel);
1313
// create an empty V matrix
1315
js.jobtype="shellcommand";
1316
js.arguments["command"]=str(format("vbmm2 -zeros %s.V %d %d")%stem%orderg%orderg);
1317
js.name="GLM-Vcreate";
1319
int n_vcreate=js.jnum;
1322
// build V matrix (V=KKt)
1323
set<int32> n_vpieces;
1328
rowfinish=rowstart+rowsperjob;
1329
if (rowfinish > orderg-1)
1332
js.jobtype="shellcommand";
1333
js.arguments["command"]=str(format("vbmm2 -xyt %1%.K %1%.K %1%.V %2% %3%")%stem%rowstart%rowfinish);
1334
sprintf(tmp,"GLM-V%d",i++);
1336
js.waitfor.insert(n_k);
1337
js.waitfor.insert(n_vcreate);
1339
n_vpieces.insert(js.jnum);
1340
n_v=js.jnum; // in case it's just one piece
1343
rowstart+=rowsperjob+1;
1344
} while (rowstart < orderg);
1349
js.jobtype="shellcommand";
1350
js.arguments["command"]=str(format("vbmm2 -assemblecols %s.V")%stem);
1351
js.name="GLM-MergeV";
1352
js.waitfor=n_vpieces;
1360
js.jobtype="shellcommand";
1361
js.arguments["command"]=str(format("vbmm2 -f3 %1%.V %1%.KG %1%.F3")%stem);
1362
js.name="GLM-makeF3";
1363
js.waitfor.insert(n_kg);
1364
js.waitfor.insert(n_v);
1369
// create an empty RV matrix
1371
js.jobtype="shellcommand";
1372
js.arguments["command"]=str(format("vbmm2 -zeros %s.RV %d %d")%stem%orderg%orderg);
1373
js.name="GLM-RVcreate";
1375
int n_rvcreate=js.jnum;
1379
set<int32> n_rvpieces;
1383
rowfinish=rowstart+rowsperjob;
1384
if (rowfinish > orderg-1)
1387
js.jobtype="shellcommand";
1388
js.arguments["command"]=str(format("vbmm2 -xyt %1%.R %1%.V %1%.RV %2% %3%")%stem%rowstart%rowfinish);
1389
sprintf(tmp,"GLM-RV%d",i);
1391
js.waitfor.insert(n_r);
1392
js.waitfor.insert(n_v);
1393
js.waitfor.insert(n_rvcreate);
1397
n_rvpieces.insert(js.jnum);
1398
n_rv=js.jnum; // in case it's one piece
1400
rowstart+=rowsperjob+1;
1401
} while (rowstart < orderg);
1406
js.jobtype="shellcommand";
1407
js.arguments["command"]=str(format("vbmm2 -assemblecols %s.RV")%stem);
1408
js.name="GLM-MergeRV";
1409
js.waitfor=n_rvpieces;
1415
// create an empty RVRV matrix
1417
js.jobtype="shellcommand";
1418
js.arguments["command"]=str(format("vbmm2 -zeros %s.RVRV %d %d")%stem%orderg%orderg);
1419
js.name="GLM-RVRVcreate";
1421
int n_rvrvcreate=js.jnum;
1425
set<int32> n_rvrvpieces;
1429
rowfinish=rowstart+rowsperjob;
1430
if (rowfinish > orderg-1)
1433
js.jobtype="shellcommand";
1434
js.arguments["command"]=str(format("vbmm2 -xy %1%.RV %1%.RV %1%.RVRV %2% %3%")%stem%rowstart%rowfinish);
1435
sprintf(tmp,"GLM-RVRV%d",i);
1437
js.waitfor.insert(n_rvrvcreate);
1438
js.waitfor.insert(n_rv);
1440
n_rvrvpieces.insert(js.jnum);
1444
rowstart+=rowsperjob+1;
1445
} while (rowstart < orderg);
1450
js.jobtype="shellcommand";
1451
js.arguments["command"]=str(format("vbmm2 -assemblecols %s.RVRV")%stem);
1452
js.name="GLM-MergeRVRV";
1453
js.waitfor=n_rvrvpieces;
1461
js.jobtype="shellcommand";
1462
js.arguments["command"]=str(format("comptraces -m %s")%stem);
1463
js.name="GLM-traces";
1464
js.waitfor.insert(n_rvrv);
1466
int n_traces=js.jnum;
1469
// tes regression steps
1470
set<int32> n_regstep;
1474
if (this->driftcorrect)
1476
for (i=0; i<pieces; i++) {
1478
js.jobtype="shellcommand";
1479
js.arguments["command"]=str(format("vbregress %s -p %d %d %s")
1480
%stem%(i+1)%pieces%flags);
1481
sprintf(tmp,"regress(%d/%d)",i+1,pieces);
1483
js.waitfor.insert(n_traces);
1484
js.waitfor.insert(n_f1);
1485
js.waitfor.insert(n_r);
1486
js.waitfor.insert(n_exofilt);
1488
n_regstep.insert(js.jnum);
1495
js.jobtype="shellcommand";
1496
js.arguments["command"]=str(format("vbmerge4d %1%.prm_part_* -o %1%.prm")%stem);
1497
js.name="mergeparams";
1498
js.waitfor=n_regstep;
1504
js.jobtype="shellcommand";
1505
js.arguments["command"]=str(format("rm %s.prm_part_*")%stem);
1506
js.name="mergeparams";
1507
js.waitfor.insert(jobnum-1);
1513
js.jobtype="shellcommand";
1514
js.arguments["command"]=str(format("vbmerge4d %1%.res_part_* -o %1%.res")%stem);
1515
js.name="mergeparams";
1516
js.waitfor=n_regstep;
1522
js.jobtype="shellcommand";
1523
js.arguments["command"]=str(format("rm %s.res_part_*")%stem);
1524
js.name="mergeparams";
1525
js.waitfor.insert(jobnum-1);
1530
// make these four jobs stand in for regstep in the waitfor chain
1532
n_regstep.insert(jobnum-4);
1533
n_regstep.insert(jobnum-3);
1534
n_regstep.insert(jobnum-2);
1535
n_regstep.insert(jobnum-1);
1538
// create smoothness estimate
1540
js.jobtype="shellcommand";
1541
js.arguments["command"]=str(format("vbse %1%.res %1%.se")%stem);
1543
js.waitfor=n_regstep;
1553
js.jobtype="shellcommand";
1554
js.arguments["command"]=str(format("glminfo -r %1% > %1%/audit.txt")%xdirname(stem));
1555
js.name="GLM-audit";
1556
js.waitfor=n_regstep;
1557
js.waitfor.insert(n_se);
1564
if (emailflag && email.size()) {
1566
js.jobtype="notify";
1567
js.arguments["email"]=email;
1568
js.arguments["msg"]="Your GLM has been solved.";
1569
js.magnitude=0; // set it to something!
1571
js.waitfor.insert(n_se);
1573
js.waitfor.insert(n_audit);
1574
js.waitfor.insert(n_f3);
1578
for (SMI j=seq.specmap.begin(); j!=seq.specmap.end(); j++)
1579
j->second.dirname=dirname;
1583
GLMParams::CreateGLMScript()
1586
int rowstart,rowfinish,mergeflag;
1590
// set pieces heuristically if the user didn't
1592
// FIXME is this any good?
1593
int cellcount=600000;
1594
pieces=(int)ceil((double)orderg*orderg/cellcount);
1595
if (pieces > orderg) pieces=orderg;
1596
if (pieces<1) pieces=1;
1598
int rowsperjob=orderg/pieces;
1599
// are we going to need the merge jobs?
1600
if (rowsperjob < orderg)
1605
vector<string> commandlist;
1608
tmps=(format("vbmakefilter -e %s.ExoFilt -lf %d -hf %d %s %s -t %d %f")
1609
%stem%lows%highs%(middles.size()?"-mf "+middles:"")
1610
%(kernelname.size()?"-k "+kernelname:"")
1612
commandlist.push_back(tmps);
1613
// make the noisemodel
1614
tmps=(format("vbmakefilter -i %s.IntrinCor %s -t %d %f")
1615
%stem%(noisemodel.size()?"-n "+noisemodel:"")%orderg%TR).str();
1616
commandlist.push_back(tmps);
1617
// make the KG matrix
1618
tmps=(format("makematkg -m %s")%stem).str();
1619
commandlist.push_back(tmps);
1620
// make the F1 matrix
1621
tmps=(format("vbmm2 -pinv %1%.KG %1%.F1")%stem).str();
1622
commandlist.push_back(tmps);
1623
// make the R matrix = I-(KG)(F1)
1624
tmps=(format("vbmm2 -imxy %1%.KG %1%.F1 %1%.R")%stem).str();
1625
commandlist.push_back(tmps);
1626
// create a K matrix
1627
tmps=(format("makematk -m %s")%stem).str();
1628
commandlist.push_back(tmps);
1629
// create an empty V matrix
1630
tmps=(format("vbmm2 -zeros %s.V %d %d")%stem%orderg%orderg).str();
1631
commandlist.push_back(tmps);
1633
// build V matrix (V=KKt)
1634
set<int32> n_vpieces;
1638
rowfinish=rowstart+rowsperjob;
1639
if (rowfinish > orderg-1)
1642
tmps=(format("vbmm2 -xyt %1%.K %1%.K %1%.V %2% %3%")%stem%rowstart%rowfinish).str();
1643
commandlist.push_back(tmps);
1644
rowstart+=rowsperjob+1;
1645
} while (rowstart < orderg);
1649
tmps=(format("vbmm2 -assemblecols %s.V")%stem).str();
1650
commandlist.push_back(tmps);
1654
tmps=(format("vbmm2 -f3 %1%.V %1%.KG %1%.F3")%stem).str();
1655
commandlist.push_back(tmps);
1657
// create an empty RV matrix
1658
tmps=(format("vbmm2 -zeros %s.RV %d %d")%stem%orderg%orderg).str();
1659
commandlist.push_back(tmps);
1664
rowfinish=rowstart+rowsperjob;
1665
if (rowfinish > orderg-1)
1667
tmps=(format("vbmm2 -xyt %1%.R %1%.V %1%.RV %2% %3%")%stem%rowstart%rowfinish).str();
1668
commandlist.push_back(tmps);
1669
rowstart+=rowsperjob+1;
1670
} while (rowstart < orderg);
1674
tmps=(format("vbmm2 -assemblecols %s.RV")%stem).str();
1675
commandlist.push_back(tmps);
1678
// create an empty RVRV matrix
1679
tmps=(format("vbmm2 -zeros %s.RVRV %d %d")%stem%orderg%orderg).str();
1680
commandlist.push_back(tmps);
1685
rowfinish=rowstart+rowsperjob;
1686
if (rowfinish > orderg-1)
1688
tmps=(format("vbmm2 -xy %1%.RV %1%.RV %1%.RVRV %2% %3%")%stem%rowstart%rowfinish).str();
1689
commandlist.push_back(tmps);
1690
rowstart+=rowsperjob+1;
1691
} while (rowstart < orderg);
1695
tmps=(format("vbmm2 -assemblecols %s.RVRV")%stem).str();
1696
commandlist.push_back(tmps);
1700
tmps=(format("comptraces -m %s")%stem).str();
1701
commandlist.push_back(tmps);
1703
// tes regression steps
1707
if (this->driftcorrect)
1709
for (i=0; i<pieces; i++) {
1710
tmps=(format("vbregress %s -p %d %d %s")
1711
%stem%(i+1)%pieces%flags).str();
1712
commandlist.push_back(tmps);
1716
tmps=str((format("vbmerge4d %1%.prm_part_* -o %1%.prm")%stem));
1717
commandlist.push_back(tmps);
1718
tmps=str((format("rm %s.prm_part_*")%stem));
1719
commandlist.push_back(tmps);
1721
tmps=str((format("vbmerge4d %1%.res_part_* -o %1%.res")%stem));
1722
commandlist.push_back(tmps);
1723
tmps=str((format("rm %s.res_part_*")%stem));
1724
commandlist.push_back(tmps);
1727
// create smoothness estimate
1728
tmps=str((format("vbse %1%.res %1%.se")%stem));
1729
commandlist.push_back(tmps);
1733
tmps=str((format("glminfo -r %1% > %1%/audit.txt")%xdirname(stem)));
1734
commandlist.push_back(tmps);
1740
GLMParams::WriteAndSubmitJobs()
1742
// copy a few things to the sequence
1745
for (SMI j=seq.specmap.begin(); j!=seq.specmap.end(); j++)
1746
j->second.dirname=dirname;
1747
return seq.Submit(vbp);
1751
GLMParams::WriteGLMFile(string fname)
1756
fp=fopen(fname.c_str(),"w");
1758
printf("[E] vbmakeglm: couldn't create glm file %s\n",fname.c_str());
1761
fprintf(fp,"lows %d\n",lows);
1762
fprintf(fp,"highs %d\n",highs);
1764
fprintf(fp,"middles %s\n",middles.c_str());
1765
fprintf(fp,"orderg %d\n",orderg);
1766
fprintf(fp,"pieces %d\n",pieces);
1767
fprintf(fp,"kernel %s\n",kernelname.c_str());
1768
fprintf(fp,"noisemodel %s\n",noisemodel.c_str());
1770
fprintf(fp,"makerandfxg\n");
1772
fprintf(fp,"gmatrix %s\n",gmatrix.c_str());
1774
fprintf(fp,"refname %s\n",refname.c_str());
1775
fprintf(fp,"pri %d\n",pri);
1776
fprintf(fp,"audit %s\n",(auditflag ? "yes" : "no"));
1777
fprintf(fp,"meannorm %s\n",(meannorm ? "yes" : "no"));
1778
fprintf(fp,"driftcorrect %s\n",(driftcorrect ? "yes" : "no"));
1779
fprintf(fp,"email %s\n",email.c_str());
1781
fprintf(fp,"glm %s\n",name.c_str());
1782
fprintf(fp,"dirname %s\n",dirname.c_str());
1783
for (int i=0; i<(int)scanlist.size(); i++)
1784
fprintf(fp,"scan %s\n",scanlist[i].c_str());
1785
fprintf(fp,"end\n");
1791
GLMParams::FixRelativePaths()
1793
string path=xgetcwd()+"/";
1794
dirname=xabsolutepath(dirname);
1795
kernelname=xabsolutepath(kernelname);
1796
noisemodel=xabsolutepath(noisemodel);
1797
refname=xabsolutepath(refname);
1798
gmatrix=xabsolutepath(gmatrix);
1799
for (size_t i=0; i<scanlist.size(); i++)
1800
scanlist[i]=xabsolutepath(scanlist[i]);
1806
printf("\nVoxBo vbmakeglm (v%s)\n",vbversion.c_str());
1807
printf("usage: vbmakeglm [flags] file [file ...]\n");
1809
printf(" -i ignore warnings, queue it anyway\n");
1810
printf(" -v don't queue, just validate\n");
1811
printf(" -s setup only (create filters)\n");
1812
printf(" -tr <rate> set data TR (ms)\n");
1813
printf(" -kr <rate> set HRF kernel sampling rate (ms)\n");
1814
printf(" -h i'm lazy, show me a sample .glm file\n");
1815
printf(" -n <name> only actually submit glms matching name\n");
1816
printf(" -x <filename> i'm really lazy, make me a sample .glm file\n");
1817
printf(" --run=<n> don't queue, run now on n cores\n");
1819
printf(" If --run is specified without =n, the default is the total number of\n");
1820
printf(" available cores minus 1 (or 1, if that would be 0).\n");
1825
vbmakeglm_make_sample(const string &fname)
1827
printf("Creating sample glm file in %s...",fname.c_str()); fflush(stdout);
1828
FILE *fp=fopen(fname.c_str(),"a");
1830
vbmakeglm_sample(fp);
1835
printf("failed.\n");
1840
vbmakeglm_show_sample()
1842
printf("Here's what a valid .glm file looks like, more or less:\n\n");
1844
vbmakeglm_sample(stdout);
1848
vbmakeglm_sample(FILE *fp)
1850
fprintf(fp,"###############################################################\n");
1851
fprintf(fp,"# Sample .glm file for VoxBo vbmakeglm.\n");
1852
fprintf(fp,"# At the top, put stuff that's common to all GLMS in this file.\n");
1853
fprintf(fp,"###############################################################\n");
1856
fprintf(fp,"# how many low frequencies do you want to filter out?\n");
1857
fprintf(fp,"lows 1\n\n");
1858
fprintf(fp,"# how many high frequencies do you want to filter out?\n");
1859
fprintf(fp,"highs 1\n\n");
1860
fprintf(fp,"# which random frequencies fo you want to filter out? (separate by spaces)\n");
1861
fprintf(fp,"# middles 15 16 21\n\n");
1862
fprintf(fp,"# how many total data points per voxel (typically time points)?\n");
1863
fprintf(fp,"orderg 160\n\n");
1864
fprintf(fp,"# how many pieces do you want to break the matrix operations into?\n");
1865
fprintf(fp,"pieces 10\n\n");
1866
fprintf(fp,"# specify a kernel for exogenous smoothing by convolution\n");
1867
fprintf(fp,"kernel /usr/local/VoxBo/elements/filters/Eigen1.ref\n\n");
1868
fprintf(fp,"# specify a model of intrinsic noise\n");
1869
fprintf(fp,"noisemodel /usr/local/VoxBo/elements/noisemodels/smooth_params.ref\n\n");
1870
fprintf(fp,"# specify the location of your G matrix, or just \"makerandfxg\" for makerandfxg\n");
1871
fprintf(fp,"gmatrix /data/Models/mytask.G\n");
1872
fprintf(fp,"# makerandfxg\n\n");
1873
fprintf(fp,"# the location of a reference function to copy into the GLM directory\n");
1874
fprintf(fp,"refname /data/Models/motor.ref\n\n");
1875
fprintf(fp,"# priority (1 for overnight, 2 for low, 3 for normal, 4 and 5 for various emergencies)\n");
1876
fprintf(fp,"pri 3\n\n");
1877
fprintf(fp,"# do some summary statistics on your GLM?\n");
1878
fprintf(fp,"audit yes\n\n");
1879
fprintf(fp,"# mean normalize? usually say yes only if you have multiple runs of BOLD data\n");
1880
fprintf(fp,"meannorm yes\n\n");
1881
fprintf(fp,"# correct linear drift? \n");
1882
fprintf(fp,"driftcorrect yes\n\n");
1883
fprintf(fp,"# your email address here\n");
1884
fprintf(fp,"email nobody@nowhere.com\n\n");
1887
fprintf(fp,"###############################################################\n");
1888
fprintf(fp,"# Then create paragraphs of specific stuff for each GLM.\n");
1889
fprintf(fp,"# (you can also override globals here)\n");
1890
fprintf(fp,"###############################################################\n");
1893
fprintf(fp,"glm larry-glm1\n");
1894
fprintf(fp,"dirname /data/study/larry/glm1\n");
1895
fprintf(fp,"scan /data/study/larry/larry01/larry01.tes\n");
1896
fprintf(fp,"scan /data/study/larry/larry02/larry02.tes\n");
1897
fprintf(fp,"end\n");
1900
fprintf(fp,"glm moe-glm1\n");
1901
fprintf(fp,"dirname /data/study/moe/glm1\n");
1902
fprintf(fp,"scan /data/study/moe/moe01/moe01.tes\n");
1903
fprintf(fp,"scan /data/study/moe/moe02/moe02.tes\n");
1904
fprintf(fp,"end\n");
1907
fprintf(fp,"glm shemp-glm1\n");
1908
fprintf(fp,"dirname /data/study/shemp/glm1\n");
1909
fprintf(fp,"scan /data/study/shemp/shemp01/shemp01.tes\n");
1910
fprintf(fp,"scan /data/study/shemp/shemp02/shemp02.tes\n");
1911
fprintf(fp,"end\n");