diff --git a/fsl_glm.cc b/fsl_glm.cc
index 41091065f07c1004ca4f45646baab37bd91523cb..a385542d72c4be3a9c728f3b979205ef09e86f0b 100644
--- a/fsl_glm.cc
+++ b/fsl_glm.cc
@@ -1,8 +1,8 @@
/* fsl_glm -
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 2006-2008 University of Oxford */
+ Copyright (C) 2006-2013 University of Oxford */
/* CCOPYRIGHT */
@@ -23,7 +23,7 @@ using namespace std;
// printed out as the help or usage message
string title=string("fsl_glm")+
- string("\nCopyright(c) 2004-2009, University of Oxford (Christian F. Beckmann)\n")+
+ string("\nAuthor: Christian F. Beckmann \nCopyright(C) 2006-2013 University of Oxford \n")+
string(" \n Simple GLM allowing temporal or spatial regression on either text data or images\n");
string examples="fsl_glm -i <input> -d <design> -o <output> [options]";
diff --git a/fsl_mvlm.cc b/fsl_mvlm.cc
index 096f74b9af1aaf270da10d4964c643123fcf7d4e..008fa22f32f453240523072af43fbb23f25833d2 100644
--- a/fsl_mvlm.cc
+++ b/fsl_mvlm.cc
@@ -1,70 +1,11 @@
/* fsl_glm -
- Christian F. Beckmann, FMRIB Image Analysis Group
-
- Copyright (C) 2006-2008 University of Oxford */
-
-/* Part of FSL - FMRIB's Software Library
- http://www.fmrib.ox.ac.uk/fsl
- fsl@fmrib.ox.ac.uk
-
- Developed at FMRIB (Oxford Centre for Functional Magnetic Resonance
- Imaging of the Brain), Department of Clinical Neurology, Oxford
- University, Oxford, UK
-
-
- LICENCE
-
- FMRIB Software Library, Release 4.0 (c) 2007, The University of
- Oxford (the "Software")
-
- The Software remains the property of the University of Oxford ("the
- University").
-
- The Software is distributed "AS IS" under this Licence solely for
- non-commercial use in the hope that it will be useful, but in order
- that the University as a charitable foundation protects its assets for
- the benefit of its educational and research purposes, the University
- makes clear that no condition is made or to be implied, nor is any
- warranty given or to be implied, as to the accuracy of the Software,
- or that it will be suitable for any particular purpose or for use
- under any specific conditions. Furthermore, the University disclaims
- all responsibility for the use which is made of the Software. It
- further disclaims any liability for the outcomes arising from using
- the Software.
-
- The Licensee agrees to indemnify the University and hold the
- University harmless from and against any and all claims, damages and
- liabilities asserted by third parties (including claims for
- negligence) which arise directly or indirectly from the use of the
- Software or the sale of any products based on the Software.
-
- No part of the Software may be reproduced, modified, transmitted or
- transferred in any form or by any means, electronic or mechanical,
- without the express permission of the University. The permission of
- the University is not required if the said reproduction, modification,
- transmission or transference is done without financial return, the
- conditions of this Licence are imposed upon the receiver of the
- product, and all original and amended source code is included in any
- transmitted product. You may be held legally responsible for any
- copyright infringement that is caused or encouraged by your failure to
- abide by these terms and conditions.
-
- You are not permitted under this Licence to use this Software
- commercially. Use for which any financial return is received shall be
- defined as commercial use, and includes (1) integration of all or part
- of the source code or the Software into a product for sale or license
- by or on behalf of Licensee to third parties or (2) use of the
- Software or any derivative of it for research with the final aim of
- developing software products for sale or license to a third party or
- (3) use of the Software or any derivative of it for research with the
- final aim of developing non-software products for sale or license to a
- third party, or (4) use of the Software to provide any service to an
- external organisation for which payment is received. If you are
- interested in using the Software commercially, please contact Isis
- Innovation Limited ("Isis"), the technology transfer company of the
- University, to negotiate a licence. Contact details are:
- innovation@isis.ox.ac.uk quoting reference DE/1112. */
+ Christian F. Beckmann, FMRIB Analysis Group
+
+ Copyright (C) 2006-2013 University of Oxford */
+
+/* CCOPYRIGHT */
+
//Header & includes
#include "libvis/miscplot.h"
@@ -109,7 +50,7 @@ using namespace std;
// printed out as the help or usage message
string title=string("fsl_mvlm (Version 1.0)")+
- string("\nCopyright(c) 2008, University of Oxford (Christian F. Beckmann)\n")+
+ string("\nAuthor: Christian F. Beckmann \nCopyright(C) 2006-2013 University of Oxford\n")+
string(" \n Multivariate Linear Model regression on\n")+
string(" time courses and/or 3D/4D imges using SVD (PCA), PLS, normalised PLS, \n")+
string(" CVA, SVD-CVA or MLM\n\n");
diff --git a/fsl_regfilt.cc b/fsl_regfilt.cc
index ffe19ebb5b57a314b77e8adbef2890a8bf925659..6db09b8cb896f735f76bd355e99ec9f09e1f1dea 100644
--- a/fsl_regfilt.cc
+++ b/fsl_regfilt.cc
@@ -1,10 +1,11 @@
/* fsl_regfilt -
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 2006-2011 University of Oxford / Christian F. Beckmann */
+ Copyright (C) 2006-2013 University of Oxford */
/* CCOPYRIGHT */
+
#include "libvis/miscplot.h"
#include "miscmaths/miscmaths.h"
#include "miscmaths/miscprob.h"
@@ -22,7 +23,7 @@ using namespace std;
// printed out as the help or usage message
string title=string("fsl_regfilt")+
- string("\n\n Copyright(c) 2011, University of Oxford (Christian F. Beckmann)\n")+
+ string("\nAuthor: Christian F. Beckmann \n Copyright(C) 2016-2013 University of Oxford\n")+
string(" Data de-noising by regressing out part of a design matrix\n")+
string(" using simple OLS regression on 4D images");
string examples="fsl_regfilt -i <input> -d <design> -f <component numbers or filter threshold> -o <out> [options]";
diff --git a/fsl_sbca.cc b/fsl_sbca.cc
index 12c1a19aee028870bf6fca6523a3c5c583e5c9b4..275ff30db25daa39f27f1acde73c84ba0efe59a8 100644
--- a/fsl_sbca.cc
+++ b/fsl_sbca.cc
@@ -1,70 +1,10 @@
/* fsl_glm -
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 2008 University of Oxford */
+ Copyright (C) 2008-2013 University of Oxford */
-/* Part of FSL - FMRIB's Software Library
- http://www.fmrib.ox.ac.uk/fsl
- fsl@fmrib.ox.ac.uk
-
- Developed at FMRIB (Oxford Centre for Functional Magnetic Resonance
- Imaging of the Brain), Department of Clinical Neurology, Oxford
- University, Oxford, UK
-
-
- LICENCE
-
- FMRIB Software Library, Release 4.0 (c) 2007, The University of
- Oxford (the "Software")
-
- The Software remains the property of the University of Oxford ("the
- University").
-
- The Software is distributed "AS IS" under this Licence solely for
- non-commercial use in the hope that it will be useful, but in order
- that the University as a charitable foundation protects its assets for
- the benefit of its educational and research purposes, the University
- makes clear that no condition is made or to be implied, nor is any
- warranty given or to be implied, as to the accuracy of the Software,
- or that it will be suitable for any particular purpose or for use
- under any specific conditions. Furthermore, the University disclaims
- all responsibility for the use which is made of the Software. It
- further disclaims any liability for the outcomes arising from using
- the Software.
-
- The Licensee agrees to indemnify the University and hold the
- University harmless from and against any and all claims, damages and
- liabilities asserted by third parties (including claims for
- negligence) which arise directly or indirectly from the use of the
- Software or the sale of any products based on the Software.
-
- No part of the Software may be reproduced, modified, transmitted or
- transferred in any form or by any means, electronic or mechanical,
- without the express permission of the University. The permission of
- the University is not required if the said reproduction, modification,
- transmission or transference is done without financial return, the
- conditions of this Licence are imposed upon the receiver of the
- product, and all original and amended source code is included in any
- transmitted product. You may be held legally responsible for any
- copyright infringement that is caused or encouraged by your failure to
- abide by these terms and conditions.
-
- You are not permitted under this Licence to use this Software
- commercially. Use for which any financial return is received shall be
- defined as commercial use, and includes (1) integration of all or part
- of the source code or the Software into a product for sale or license
- by or on behalf of Licensee to third parties or (2) use of the
- Software or any derivative of it for research with the final aim of
- developing software products for sale or license to a third party or
- (3) use of the Software or any derivative of it for research with the
- final aim of developing non-software products for sale or license to a
- third party, or (4) use of the Software to provide any service to an
- external organisation for which payment is received. If you are
- interested in using the Software commercially, please contact Isis
- Innovation Limited ("Isis"), the technology transfer company of the
- University, to negotiate a licence. Contact details are:
- innovation@isis.ox.ac.uk quoting reference DE/1112. */
+/* CCOPYRIGHT */
#include "libvis/miscplot.h"
#include "miscmaths/miscmaths.h"
@@ -83,7 +23,7 @@ using namespace std;
// printed out as the help or usage message
string title=string("fsl_sbca (Version 1.0)")+
- string("\nCopyright(c) 2008, University of Oxford (Christian F. Beckmann)\n")+
+ string("\nAuthor: Christian F. Beckmann \nCopyright(C) 2008-2013 University of Oxford \n")+
string(" \n Performs seed-based correlation analysis on FMRI data\n")+
string(" using either a single seed coordinate or a seed mask \n")+
string(" ");
diff --git a/fsl_schurprod.cc b/fsl_schurprod.cc
index 8b46bad5e01313625427ddc7795707a4b0e48478..04bfc7b8c2a90ef9b50182eae03feee625a8971e 100644
--- a/fsl_schurprod.cc
+++ b/fsl_schurprod.cc
@@ -1,70 +1,10 @@
/* fsl_glm -
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 2008 University of Oxford */
+ Copyright (C) 2008-2013 University of Oxford */
-/* Part of FSL - FMRIB's Software Library
- http://www.fmrib.ox.ac.uk/fsl
- fsl@fmrib.ox.ac.uk
-
- Developed at FMRIB (Oxford Centre for Functional Magnetic Resonance
- Imaging of the Brain), Department of Clinical Neurology, Oxford
- University, Oxford, UK
-
-
- LICENCE
-
- FMRIB Software Library, Release 4.0 (c) 2007, The University of
- Oxford (the "Software")
-
- The Software remains the property of the University of Oxford ("the
- University").
-
- The Software is distributed "AS IS" under this Licence solely for
- non-commercial use in the hope that it will be useful, but in order
- that the University as a charitable foundation protects its assets for
- the benefit of its educational and research purposes, the University
- makes clear that no condition is made or to be implied, nor is any
- warranty given or to be implied, as to the accuracy of the Software,
- or that it will be suitable for any particular purpose or for use
- under any specific conditions. Furthermore, the University disclaims
- all responsibility for the use which is made of the Software. It
- further disclaims any liability for the outcomes arising from using
- the Software.
-
- The Licensee agrees to indemnify the University and hold the
- University harmless from and against any and all claims, damages and
- liabilities asserted by third parties (including claims for
- negligence) which arise directly or indirectly from the use of the
- Software or the sale of any products based on the Software.
-
- No part of the Software may be reproduced, modified, transmitted or
- transferred in any form or by any means, electronic or mechanical,
- without the express permission of the University. The permission of
- the University is not required if the said reproduction, modification,
- transmission or transference is done without financial return, the
- conditions of this Licence are imposed upon the receiver of the
- product, and all original and amended source code is included in any
- transmitted product. You may be held legally responsible for any
- copyright infringement that is caused or encouraged by your failure to
- abide by these terms and conditions.
-
- You are not permitted under this Licence to use this Software
- commercially. Use for which any financial return is received shall be
- defined as commercial use, and includes (1) integration of all or part
- of the source code or the Software into a product for sale or license
- by or on behalf of Licensee to third parties or (2) use of the
- Software or any derivative of it for research with the final aim of
- developing software products for sale or license to a third party or
- (3) use of the Software or any derivative of it for research with the
- final aim of developing non-software products for sale or license to a
- third party, or (4) use of the Software to provide any service to an
- external organisation for which payment is received. If you are
- interested in using the Software commercially, please contact Isis
- Innovation Limited ("Isis"), the technology transfer company of the
- University, to negotiate a licence. Contact details are:
- innovation@isis.ox.ac.uk quoting reference DE/1112. */
+/* CCOPYRIGHT */
#include "libvis/miscplot.h"
#include "miscmaths/miscmaths.h"
@@ -83,7 +23,7 @@ using namespace std;
// printed out as the help or usage message
string title=string("fsl_schurprod (Version 1.0)")+
- string("\nCopyright(c) 2011, Christian F. Beckmann\n")+
+ string("\nAuthor: Christian F. Beckmann \nCopyright(C) University of Oxford\n")+
string(" \n Generates element-wise matrix products or product of matrices against vectors from 4D data\n")+
string(" ");
string examples="fsl_schurprod -i <input> -d <time series> -o <basename> [options]";
diff --git a/ggmix.cc b/ggmix.cc
index 4b77a6ff4866d1d09b852f8821d46500d4540a77..9e2cd7f8c9934c614764f0857ca6805c7dfaa95e 100644
--- a/ggmix.cc
+++ b/ggmix.cc
@@ -3,9 +3,9 @@
ggmix.cc - Gaussian & Gaussian/Gamma Mixture Model
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-20013 University of Oxford */
/* CCOPYRIGHT */
diff --git a/ggmix.h b/ggmix.h
index e97a4e157beaa2b72471399efc259f72bb09db2a..716c6902e484e6828cdcb3b57700f500512f41b2 100644
--- a/ggmix.h
+++ b/ggmix.h
@@ -3,9 +3,9 @@
ggmix.h - class for Gaussian/Gamma Mixture Model
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
diff --git a/meldata.cc b/meldata.cc
index 23377a5b16d5af169513d529c2f90103db5511f9..18bdc4c42ca582bf0db12af13f4cbc646b05ed6e 100644
--- a/meldata.cc
+++ b/meldata.cc
@@ -27,6 +27,8 @@ namespace Melodic{
ReturnMatrix MelodicData::process_file(string fname, int numfiles)
{
+ dbgmsg(string("START: process_file") << endl);
+
Matrix tmpData;
{
volume4D<float> RawData;
@@ -97,21 +99,24 @@ namespace Melodic{
}
//variance - normalisation
- memmsg(" before VN ");
+
if(opts.varnorm.value()){
+ memmsg(" before VN ");
message(" Normalising by voxel-wise variance ...");
- if(stdDev.Storage()==0)
+ outMsize("stdDev",stdDev);
+// if(stdDev.Storage()==0)
stdDev = varnorm(tmpData,std::min(30,tmpData.Nrows()-1),
- opts.vn_level.value())/numfiles;
- else
- stdDev += varnorm(tmpData,std::min(30,tmpData.Nrows()-1),
- opts.vn_level.value())/numfiles;
+ opts.vn_level.value(), opts.econ.value());
+// else
+// stdDev += varnorm(tmpData,std::min(30,tmpData.Nrows()-1),
+// opts.vn_level.value(), opts.econ.value())/numfiles;
stdDevi = pow(stdDev,-1);
memmsg(" in VN ");
message(" done" << endl);
}
tmpData.Release();
+ dbgmsg(string("END: process_file") << endl);
return tmpData;
}
@@ -240,11 +245,12 @@ namespace Melodic{
basicGLM tmpglm;
for(int ctr = 0; ctr < numfiles; ctr++){
tmpData = process_file(opts.inputfname.value().at(ctr), numfiles);
- outMsize("IC",IC);
//may want to remove the spatial means first
tmpglm.olsfit(remmean(tmpData.t(),1),remmean(IC.t(),1),tmpcont);
- s1=tmpglm.get_beta();
+ s1=tmpglm.get_beta().t();
+ outMsize("s1",s1);
+ outMsize("alltcs",alltcs);
if(alltcs.Storage()==0)
alltcs=s1;
else
@@ -252,12 +258,16 @@ namespace Melodic{
// output DR
if(opts.dr_out.value()){
- write_ascii_matrix(drO.appendDir("dr_stage1_subject"+num2str(ctr)+".txt"),s1);
+
+ dbgmsg(string("START: dual_regression output") << endl);
+ write_ascii_matrix(drO.appendDir("dr_stage1_subject"+num2str(ctr,4)+".txt"),s1);
//des_norm
s1 = SP(s1,ones(s1.Nrows(),1)*pow(stdev(s1,1),-1));
tmpglm.olsfit(remmean(tmpData),remmean(s1,1),tmpcont);
+ s2=tmpglm.get_beta();
+ save4D(s2,string("dr/dr_stage2_subject"+num2str(ctr,4)));
s2=tmpglm.get_z();
- save4D(s2,string("dr/dr_stage2_subject"+num2str(ctr)));
+ save4D(s2,string("dr/dr_stage2_subject"+num2str(ctr,4)+"_Z"));
}
}
@@ -265,6 +275,21 @@ namespace Melodic{
tmpcont << alltcs.Column(ctr);
add_Tmodes(tmpcont);
}
+
+ for(int ctrC = 1; ctrC <=IC.Nrows(); ctrC++){
+ Matrix tmpall = zeros(numfiles,IC.Ncols());
+ string fnout = string("dr/dr_stage2_ic"+num2str(ctrC-1,4));
+ for(int ctrS = 0; ctrS < numfiles; ctrS++){
+ string fnin = logger.appendDir(string("dr/dr_stage2_subject"+num2str(ctrS,4)));
+ dbgmsg(fnout << endl << fnin << endl);
+ volume4D<float> vol;
+ read_volume4DROI(vol,fnin,0,0,0,ctrC-1,-1,-1,-1,ctrC-1);
+
+ Matrix tmp2 = vol.matrix(Mask);
+ tmpall.Row(ctrS+1) << vol.matrix(Mask);
+ }
+ save4D(tmpall,fnout);
+ }
opts.varnorm.set_T(tmpvarnorm);
dbgmsg(string("END: dual_regression") << endl);
@@ -283,7 +308,7 @@ namespace Melodic{
void MelodicData::setup_classic()
{
-
+ dbgmsg(string("START: setup_classic") << endl);
Matrix alldat, tmpData;
bool tmpvarnorm = opts.varnorm.value();
@@ -298,14 +323,10 @@ namespace Melodic{
exit(2);
}
- if(opts.debug.value())
- save4D(alldat,string("preproc_dat") + num2str(1));
for(int ctr = 1; ctr < numfiles; ctr++){
tmpData = process_file(opts.inputfname.value().at(ctr), numfiles) / numfiles;
- if(opts.debug.value())
- save4D(tmpData,string("preproc_dat") + num2str(ctr+1));
if(tmpData.Ncols() == alldat.Ncols() && tmpData.Nrows() == alldat.Nrows())
- alldat += tmpData;
+ alldat = alldat + tmpData;
else{
if(opts.approach.value()==string("tica")){
cerr << "ERROR:: data dimensions do not match, TICA not possible \n\n";
@@ -333,7 +354,7 @@ namespace Melodic{
if((numfiles > 1 ) && opts.joined_vn.value() && tmpvarnorm){
//variance - normalisation
message(endl<<"Normalising jointly by voxel-wise variance ...");
- stdDev = varnorm(alldat,alldat.Nrows(),3.1);
+ stdDev = varnorm(alldat,alldat.Nrows(),opts.vn_level.value(),opts.econ.value());
stdDevi = pow(stdDev,-1);
message(" done" << endl);
}
@@ -395,7 +416,7 @@ namespace Melodic{
else {
dbgmsg("Multi-Subject ICA");
- stdDev.CleanUp();
+ //stdDev.CleanUp();
for(int ctr = 0; ctr < numfiles; ctr++){
tmpData = process_file(opts.inputfname.value().at(ctr), numfiles);
@@ -409,11 +430,11 @@ namespace Melodic{
Matrix newWM,newDWM;
if(!opts.joined_whiten.value()){
message(" Individual whitening in a " << order << " dimensional subspace " << endl);
- std_pca(tmpData, RXweight, Corr, pcaE, pcaD);
+ std_pca(tmpData, RXweight, Corr, pcaE, pcaD, opts.econ.value());
calc_white(pcaE, pcaD, order, newWM, newDWM);
}else{
if(!opts.dr_pca.value()){
- std_pca(whiteMatrix*tmpData, RXweight, Corr, pcaE, pcaD);
+ std_pca(whiteMatrix*tmpData, RXweight, Corr, pcaE, pcaD, opts.econ.value());
calc_white(pcaE, pcaD, order, newWM, newDWM);
newDWM=(dewhiteMatrix*newDWM);
newWM=(newWM*whiteMatrix);
@@ -426,7 +447,7 @@ namespace Melodic{
remmean(tmp1,2);
tmp1 *= tmpData.t();
tmp2 = pinv(tmp1.t()).t();
- std_pca(tmp1 * tmpData, RXweight, Corr, pcaE, pcaD);
+ std_pca(tmp1 * tmpData, RXweight, Corr, pcaE, pcaD, opts.econ.value());
calc_white(pcaE, pcaD, order, newWM, newDWM);
newDWM=(tmp2*newDWM);
newWM=(newWM * tmp1);
@@ -444,12 +465,13 @@ namespace Melodic{
}
}
opts.varnorm.set_T(tmpvarnorm);
-
+ dbgmsg(string("END: setup_classic") << endl);
+
}
void MelodicData::setup_migp()
{
- dbgmsg("starting MIGP");
+ dbgmsg(string("START: setup_migp") << endl);
std::vector<int> myctr;
for (int i=0; i< numfiles ; ++i) myctr.push_back(i);
@@ -486,7 +508,7 @@ namespace Melodic{
message(" Reducing data matrix to a " << opt.migpN.value() << " dimensional subspace " << endl);
Matrix pcaE, Corr;
RowVector pcaD;
- std_pca(Data, RXweight, Corr, pcaE, pcaD);
+ std_pca(Data, RXweight, Corr, pcaE, pcaD, opts.econ.value());
pcaE = pcaE.Columns(pcaE.Ncols()-opts.migpN.value()+1,pcaE.Ncols());
Data = pcaE.t() * Data;
}
@@ -506,10 +528,12 @@ namespace Melodic{
DWM.push_back(tmp);
WM.push_back(tmp);
opts.varnorm.set_T(tmpvarnorm);
+ dbgmsg(string("END: setup_migp") << endl);
}
void MelodicData::setup()
{
+ dbgmsg(string("START: setup") << endl);
numfiles = (int)opts.inputfname.value().size();
setup_misc();
@@ -536,10 +560,13 @@ namespace Melodic{
//save the mean & mask
save_volume(Mask,logger.appendDir("mask"));
save_volume(Mean,logger.appendDir("mean"));
+
+ dbgmsg(string("START: setup") << endl);
} // void setup()
void MelodicData::setup_misc()
{
+ dbgmsg(string("START: setup_misc") << endl);
//initialize Mean
read_volume(Mean,opts.inputfname.value().at(0));
@@ -617,6 +644,9 @@ namespace Melodic{
}
}
remmean(Tdes);
+
+ dbgmsg(string("END: setup_misc") << endl);
+
}
void MelodicData::save()
diff --git a/meldata.h b/meldata.h
index 0ea39694f6f9c381ac11c3f52883b32e6def7872..8103868d450c743e92a9f63fc856148579ab9c4f 100644
--- a/meldata.h
+++ b/meldata.h
@@ -3,9 +3,9 @@
meldata.h - data container class
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
diff --git a/melgmix.cc b/melgmix.cc
index 12e6f425dcc80acd9647328b0b3cc095fef8e89b..c8a8707570c6ee5b793fb820f59fd33383a61252 100644
--- a/melgmix.cc
+++ b/melgmix.cc
@@ -3,9 +3,9 @@
melgmix.cc - Gaussian Mixture Model
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
diff --git a/melgmix.h b/melgmix.h
index 793d3660d4df5a291db0ad611249474985260f2e..58e06e5f73e6a139490466ef8e7cf2297e229e5a 100644
--- a/melgmix.h
+++ b/melgmix.h
@@ -3,9 +3,9 @@
melgmix.h - class for Gaussian/Gamma Mixture Model
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
diff --git a/melhlprfns.cc b/melhlprfns.cc
index 3692b7d5867f6da87abe332de1c3dd0624558bad..87a1c46012575f187be441cc0eea22bf521aa131 100644
--- a/melhlprfns.cc
+++ b/melhlprfns.cc
@@ -3,9 +3,9 @@
melhlprfns.cc - misc functions
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
@@ -125,12 +125,12 @@ namespace Melodic{
return meanimg;
} //void convert_to_pbsc
- RowVector varnorm(Matrix& in, int dim, float level)
+ RowVector varnorm(Matrix& in, int dim, float level, int econ)
{
Matrix Corr;
- Corr = calc_corr(in);
+ Corr = calc_corr(in,econ);
RowVector out;
- out = varnorm(in,Corr,dim,level);
+ out = varnorm(in,Corr,dim,level, econ);
return out;
} //RowVector varnorm
@@ -141,13 +141,13 @@ namespace Melodic{
}
}
- RowVector varnorm(Matrix& in, Matrix& Corr, int dim, float level)
+ RowVector varnorm(Matrix& in, Matrix& Corr, int dim, float level, int econ)
{
Matrix tmpE, white, dewhite;
RowVector tmpD, tmpD2;
- std_pca(remmean(in,2), Corr, tmpE, tmpD);
+ std_pca(remmean(in,2), Corr, tmpE, tmpD, econ);
calc_white(tmpE,tmpD, dim, white, dewhite);
Matrix ws = white * in;
@@ -166,11 +166,11 @@ namespace Melodic{
return tmpD;
} //RowVector varnorm
- Matrix SP2(const Matrix& in, const Matrix& weights, bool econ)
+ Matrix SP2(const Matrix& in, const Matrix& weights, int econ)
{
Matrix Res;
Res = in;
- if(econ){
+ if(econ>0){
ColumnVector tmp;
for(int ctr=1; ctr <= in.Ncols(); ctr++){
tmp = in.Column(ctr);
@@ -211,16 +211,37 @@ namespace Melodic{
return out;
}
- Matrix calc_corr(const Matrix& in, bool econ)
- {
+ Matrix calc_corr(const Matrix& in, int econ)
+ {
Matrix Res;
- Res = zeros(in.Nrows(),in.Nrows());
- if(econ){
- ColumnVector tmp;
- for(int ctr=1; ctr <= in.Ncols(); ctr++){
- tmp = in.Column(ctr);
- tmp = tmp - mean(tmp).AsScalar();
- Res += (tmp * tmp.t()) / in.Ncols();
+ int nrows=in.Nrows();
+ int ncols=in.Ncols();
+ Res = zeros(nrows,nrows);
+
+ if(econ>0){
+ RowVector colmeans(ncols);
+ for (int n=1; n<=ncols; n++) {
+ colmeans(n)=0;
+ for (int m=1; m<=nrows; m++) {
+ colmeans(n)+=in(m,n);
+ }
+ colmeans(n)/=nrows;
+ }
+ int dcol = econ;
+ Matrix suba;
+
+ for(int ctr=1; ctr <= in.Ncols(); ctr+=dcol){
+ suba=in.SubMatrix(1,nrows,ctr,Min(ctr+dcol-1,ncols));
+ int scolmax = suba.Ncols();
+
+ for (int n=1; n<=scolmax; n++) {
+ double cmean=colmeans(ctr + n - 1);
+ for (int m=1; m<=nrows; m++) {
+ suba(m,n)-=cmean;
+ }
+ }
+
+ Res += suba*suba.t() / ncols;
}
}
else
@@ -228,23 +249,43 @@ namespace Melodic{
return Res;
} //Matrix calc_corr
- Matrix calc_corr(const Matrix& in, const Matrix& weights, bool econ)
+ Matrix calc_corr(const Matrix& in, const Matrix& weights, int econ)
{
Matrix Res;
- Res = zeros(in.Nrows(),in.Nrows());
- Matrix localweights;
+ int nrows=in.Nrows();
+ int ncols=in.Ncols();
+ Res = zeros(nrows,nrows);
+ RowVector localweights;
if(weights.Storage() == 0)
- localweights = ones(1,in.Ncols());
+ localweights = ones(1,ncols);
else
localweights = weights;
- if(econ){
- ColumnVector tmp;
- for(int ctr=1; ctr <= in.Ncols(); ctr++){
- tmp = in.Column(ctr);
- tmp = tmp - mean(tmp).AsScalar();
- tmp = tmp * localweights(1,ctr);
- Res += (tmp * tmp.t()) / in.Ncols();
- }
+ if(econ>0){
+ RowVector colmeans(ncols);
+ for (int n=1; n<=ncols; n++) {
+ colmeans(n)=0;
+ for (int m=1; m<=nrows; m++) {
+ colmeans(n)+=in(m,n);
+ }
+ colmeans(n)/=nrows;
+ }
+ int dcol = econ;
+ Matrix suba;
+
+ for(int ctr=1; ctr <= in.Ncols(); ctr+=dcol){
+ suba=in.SubMatrix(1,nrows,ctr,Min(ctr+dcol-1,ncols));
+ int scolmax = suba.Ncols();
+
+ for (int n=1; n<=scolmax; n++) {
+ double cmean=colmeans(ctr + n - 1);
+ double locw =localweights(ctr + n - 1);
+ for (int m=1; m<=nrows; m++) {
+ suba(m,n)-=cmean;
+ suba(m,n)*= locw;
+ }
+ }
+ Res += suba*suba.t() / ncols;
+ }
}
else{
Res = SP2(in,localweights);
@@ -253,6 +294,8 @@ namespace Melodic{
return Res;
} //Matrix calc_corr
+
+
float calc_white(const Matrix& tmpE, const RowVector& tmpD, const RowVector& PercEV, int dim, Matrix& param, Matrix& paramS, Matrix& white, Matrix& dewhite)
{
@@ -343,12 +386,12 @@ namespace Melodic{
} //Matrix calc_white
- void std_pca(const Matrix& Mat, const Matrix& weights, Matrix& Corr, Matrix& evecs, RowVector& evals)
+ void std_pca(const Matrix& Mat, const Matrix& weights, Matrix& Corr, Matrix& evecs, RowVector& evals, int econ)
{
if(weights.Storage()>0)
- Corr << calc_corr(Mat, weights);
+ Corr << calc_corr(Mat, weights, econ);
else
- Corr << calc_corr(Mat);
+ Corr << calc_corr(Mat, econ);
SymmetricMatrix tmp;
tmp << Corr;
DiagonalMatrix tmpD;
@@ -356,10 +399,10 @@ namespace Melodic{
evals = tmpD.AsRow();
} //void std_pca
- void std_pca(const Matrix& Mat, Matrix& Corr, Matrix& evecs, RowVector& evals)
+ void std_pca(const Matrix& Mat, Matrix& Corr, Matrix& evecs, RowVector& evals, int econ)
{
Matrix weights;
- std_pca(Mat,weights,Corr,evecs,evals);
+ std_pca(Mat,weights,Corr,evecs,evals, econ);
} //void std_pca
void em_pca(const Matrix& Mat, Matrix& evecs, RowVector& evals, int num_pc, int iter)
diff --git a/melhlprfns.h b/melhlprfns.h
index 5903ed2638e039cf0e5ef8c86f6ba1f4516eabb4..1d2942a3211a7d79ce8e3255d683f25c9006c8ad 100644
--- a/melhlprfns.h
+++ b/melhlprfns.h
@@ -3,9 +3,9 @@
melhlprfns.cc - misc functions
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
@@ -44,11 +44,11 @@ namespace Melodic{
Matrix convert_to_pbsc(Matrix& Mat);
- RowVector varnorm(Matrix& in, int dim = 30, float level = 1.6);
+ RowVector varnorm(Matrix& in, int dim = 30, float level = 1.6, int econ = 0);
void varnorm(Matrix& in, const RowVector& vars);
- RowVector varnorm(Matrix& in, Matrix& Corr, int dim = 30, float level = 1.6);
+ RowVector varnorm(Matrix& in, Matrix& Corr, int dim = 30, float level = 1.6, int econ = 0);
- Matrix SP2(const Matrix& in, const Matrix& weights, bool econ = 0);
+ Matrix SP2(const Matrix& in, const Matrix& weights, int econ = 20000);
void SP3(Matrix& in, const Matrix& weights);
RowVector Feta(int n1,int n2);
@@ -56,8 +56,8 @@ namespace Melodic{
Matrix corrcoef(const Matrix& in1, const Matrix& in2);
Matrix corrcoef(const Matrix& in1, const Matrix& in2, const Matrix& part);
- Matrix calc_corr(const Matrix& in, bool econ = 1);
- Matrix calc_corr(const Matrix& in, const Matrix& weights, bool econ = 1);
+ Matrix calc_corr(const Matrix& in, int econ = 20000);
+ Matrix calc_corr(const Matrix& in, const Matrix& weights, int econ = 20000);
float calc_white(const Matrix& tmpE, const RowVector& tmpD, const RowVector& PercEV, int dim, Matrix& param, Matrix& paramS, Matrix& white, Matrix& dewhite);
float calc_white(const Matrix& tmpE, const RowVector& tmpD, const RowVector& PercEV, int dim, Matrix& white, Matrix& dewhite);
@@ -65,8 +65,8 @@ namespace Melodic{
void calc_white(const Matrix& tmpE, const RowVector& tmpD, int dim, Matrix& white, Matrix& dewhite);
void calc_white(const Matrix& Corr, int dim, Matrix& white, Matrix& dewhite);
- void std_pca(const Matrix& Mat, Matrix& Corr, Matrix& evecs, RowVector& evals);
- void std_pca(const Matrix& Mat, const Matrix& weights, Matrix& Corr, Matrix& evecs, RowVector& evals);
+ void std_pca(const Matrix& Mat, Matrix& Corr, Matrix& evecs, RowVector& evals, int econ = 20000);
+ void std_pca(const Matrix& Mat, const Matrix& weights, Matrix& Corr, Matrix& evecs, RowVector& evals, int econ = 20000);
void em_pca(const Matrix& Mat, Matrix& evecs, RowVector& evals, int num_pc = 1, int iter = 20);
void em_pca(const Matrix& Mat, Matrix& guess, Matrix& evecs, RowVector& evals, int num_pc = 1, int iter = 20);
@@ -92,7 +92,7 @@ namespace Melodic{
ColumnVector gen_ar(const ColumnVector& in, int maxorder = 1);
Matrix gen_ar(const Matrix& in, int maxorder);
Matrix gen_arCorr(const Matrix& in, int maxorder);
-
+
class basicGLM{
public:
diff --git a/melica.cc b/melica.cc
index e6a14583b34c89d91060dcebe7a24e9f25bd53a0..d517e0c67ceba97f3b519f585ae8bcc02e47a626 100644
--- a/melica.cc
+++ b/melica.cc
@@ -3,9 +3,9 @@
melica.cc - ICA estimation
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
diff --git a/melica.h b/melica.h
index af098f38c3b2b5cc33dadd27edd95e2747b1ce20..dc9594033322e511a137af4d90e28eaedd890f09 100644
--- a/melica.h
+++ b/melica.h
@@ -3,9 +3,9 @@
melica.h - ICA estimation
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
diff --git a/melodic.cc b/melodic.cc
index f1a4558f6be819a885cb43fbfd32e9bfb306ed74..bf7bde542a94104cc23c8642c5b05d90095c9b17 100644
--- a/melodic.cc
+++ b/melodic.cc
@@ -3,7 +3,7 @@
melodic.cc - main program file
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
Copyright (C) 1999-2013 University of Oxford */
diff --git a/melodic.h b/melodic.h
index 7c3776c426897455ccc9d7ac044d4eebff37b32d..2e938b58122204b2adfb16c0c80a9e88ca21d68f 100644
--- a/melodic.h
+++ b/melodic.h
@@ -3,9 +3,9 @@
melodic.h - main program header
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
@@ -66,13 +66,13 @@
namespace Melodic{
-const string version = "3.12";
+const string version = "3.13";
// The two strings below specify the title and example usage that is
// printed out as the help or usage message
const string title=string("MELODIC (Version ")+version+")"+
string("\n Multivariate Exploratory Linear Optimised Decomposition into Independent Components\n")+
- string(" Copyright(c) 2001-2008, University of Oxford (Christian F. Beckmann)");
+ string("\nAuthor: Christian F. Beckmann \n Copyright(c) 2001-2013 University of Oxford");
const string usageexmpl=string(" melodic -i <filename> <options>")+
string("\n \t \t to run melodic")+
diff --git a/meloptions.cc b/meloptions.cc
index 4efdae2d407d51220fe869b1c15603cdcd046703..3b2f67ce2361e2524cc61170422cf90c36ef26ee 100644
--- a/meloptions.cc
+++ b/meloptions.cc
@@ -3,9 +3,9 @@
meloptions.cc - class for command line options
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
diff --git a/meloptions.h b/meloptions.h
index a7cd248c177c27dc547fac5feeb8768c182a170c..65e17cd4ce5a5a248ad3056dd9f875c7af681469 100644
--- a/meloptions.h
+++ b/meloptions.h
@@ -3,9 +3,9 @@
meloptions.h - class for command line options
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
@@ -134,6 +134,7 @@ class MelodicOptions {
Option<bool> temporal;
Option<float> retryfactor;
+ Option<int> econ;
void parse_command_line(int argc, char** argv, Log& logger, const string &p_version);
@@ -194,7 +195,7 @@ class MelodicOptions {
joined_whiten(string("--sep_whiten"), false,
string("switch on separate whitening"),
false, no_argument, false),
- joined_vn(string("--group_vn"), false,
+ joined_vn(string("--sep_vn"), true,
string("switch on group variance nomalisation (as opposed to separate VN)"),
false, no_argument, false),
dr_pca(string("--mod_pca"), true,
@@ -206,14 +207,14 @@ class MelodicOptions {
migpN(string("--migpN"), 0,
string("Number of internal Eigenmaps"),
false, requires_argument, false),
- migp_shuffle(string("--migp_order"), true,
+ migp_shuffle(string("--migp_shuffle"), true,
string("Randomise MIGP file order (default: TRUE)"),
false, no_argument, false),
migp_factor(string("--migp_factor"), 2,
string("Internal Factor of mem-threshold relative to number of Eigenmaps (default: 2)"),
false, requires_argument, false),
dr(string("--dr"), false,
- string("Dual Regression (default: TRUE)"),
+ string("Dual Regression (default: false)"),
false, no_argument, false),
dr_out(string("--dr_out"), false,
string("Dual Regression output for MIGP/concat ICA"),
@@ -245,10 +246,10 @@ class MelodicOptions {
tsmooth(string("--spca"), false,
string("smooth the eigenvectors prior to IC decomposition"),
false, no_argument, false),
- epsilon(string("--eps,--epsilon"), 0.00005,
+ epsilon(string("--eps"), 0.00005,
string("minimum error change"),
false, requires_argument),
- epsilonS(string("--epsS,--epsilonS"), 0.03,
+ epsilonS(string("--epsS"), 0.03,
string("minimum error change for rank-1 approximation in TICA"),
false, requires_argument),
maxNumItt(string("--maxit"), 500,
@@ -267,16 +268,16 @@ class MelodicOptions {
string("\tswitch off mixture modelling on IC maps\n "),
false, no_argument),
ICsfname(string("--ICs"), string(""),
- string("\tfilename of the IC components file for mixture modelling"),
+ string("\tinput filename of the IC components file for mixture modelling"),
false, requires_argument),
filtermix(string("--mix"), string(""),
- string("\tmixing matrix for mixture modelling / filtering"),
+ string("\tinput filename of mixing matrix for mixture modelling / filtering"),
false, requires_argument),
smodename(string("--smode"), string(""),
- string("\tmatrix of session modes for report generation"),
+ string("\tinput filename of matrix of session modes for report generation"),
false, requires_argument),
filter(string("-f,--filter"), string(""),
- string("component numbers to remove\n "),
+ string("list of component numbers to remove\n "),
false, requires_argument),
genreport(string("--report"), false,
string("generate Melodic web report"),
@@ -407,6 +408,9 @@ class MelodicOptions {
retryfactor(string("--retryfactor"), float(0.95),
string("multiplicative factor for determining new dim if estimated dim fails to converge"),
false, requires_argument, false),
+ econ(string("--econ"), 20000,
+ string("set ctrl parameter for helperfns econ mode"),
+ false, requires_argument, false),
options(title, usageexmpl)
{
try {
@@ -494,6 +498,7 @@ class MelodicOptions {
options.add(guess_remderiv);
options.add(temporal);
options.add(retryfactor);
+ options.add(econ);
}
catch(X_OptionError& e) {
options.usage();
diff --git a/melpca.cc b/melpca.cc
index 6bf883463c4b31eda2ddb2e21ebef2dfd73f2518..6f7b6cf92c9222591ce641384ec88b70e71e893b 100644
--- a/melpca.cc
+++ b/melpca.cc
@@ -3,9 +3,9 @@
melpca.cc - PCA and whitening
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
diff --git a/melpca.h b/melpca.h
index 57ff0ed2a55751a8fc5f5f2645f03a4509583dba..ee786a306141d27e89dcbed8b58f78519ea85264 100644
--- a/melpca.h
+++ b/melpca.h
@@ -3,9 +3,9 @@
melpca.h - PCA and whitening
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
diff --git a/melreport.cc b/melreport.cc
index de801aae987dec2e91b5822e125dadf045ce02d5..417f3ef78497e2d93826c2f13e58526a5f8e276e 100644
--- a/melreport.cc
+++ b/melreport.cc
@@ -3,9 +3,9 @@
melreport.cc - report generation
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
diff --git a/melreport.h b/melreport.h
index fa6d07337748faab1ca368fc753e8c29d7aed467..c457f94829328c74a55c3a62d72f39f474d76ca4 100644
--- a/melreport.h
+++ b/melreport.h
@@ -3,9 +3,9 @@
melreport.h - report generation
- Christian F. Beckmann, FMRIB Image Analysis Group
+ Christian F. Beckmann, FMRIB Analysis Group
- Copyright (C) 1999-2008 University of Oxford */
+ Copyright (C) 1999-2013 University of Oxford */
/* CCOPYRIGHT */
diff --git a/test.cc b/test.cc
index 29e422c532d6cb8adc5335d43a7593454792bb3b..73c75585b04407dd1d2fb1a0aef7c98e0dcc9a49 100644
--- a/test.cc
+++ b/test.cc
@@ -1,6 +1,8 @@
/* test.cc
-
- Copyright (C) 1999-2008 University of Oxford */
+
+ Christian F. Beckmann, FMRIB Analysis Group
+
+ Copyright (C) 1999-20013 University of Oxford */
/* CCOPYRIGHT */
@@ -42,20 +44,29 @@ using namespace MISCPLOT;
using namespace MISCMATHS;
using namespace Utilities;
using namespace std;
+using namespace Melodic;
+
+// GLOBALS
+
+time_t sttime;
+
// The two strings below specify the title and example usage that is
// printed out as the help or usage message
string title=string("fsl_BLAH")+
- string("\nCopyright(c) 2008, University of Oxford (Christian F. Beckmann)\n")+
+ string("\nAuthor: Christian F. Beckmann \nCopyright(c) 2008-2013 University of Oxford\n")+
string(" \n \n")+
string(" \n");
string examples="fsl_BLAH [options]";
//Command line Options {
- Option<string> fnin(string("-i,--in"), string(""),
+ Option<string> fnin(string("-i,--in"), string(""),
string("input file name (matrix 3D or 4D image)"),
true, requires_argument);
+ Option<string> fnmask(string("-m"), string(""),
+ string("mask file name "),
+ true, requires_argument);
Option<int> help(string("-h,--help"), 0,
string("display this help text"),
false,no_argument);
@@ -72,37 +83,83 @@ using namespace std;
// Local functions
+void tic(){
+ sttime = time(NULL);
+}
+
+void toc(){
+
+ time_t current_time;
+ current_time = time(NULL) - sttime;
+ cerr << endl << "TOC: " << string(ctime(¤t_time)) << endl;
+}
+
+Matrix calccorr(const Matrix& in, bool econ)
+ {
+ Matrix Res;
+ Res = zeros(in.Nrows(),in.Nrows());
+ if(econ){
+ ColumnVector tmp;
+ for(int ctr=1; ctr <= in.Ncols(); ctr++){
+ tmp = in.Column(ctr);
+ tmp = tmp - mean(tmp).AsScalar();
+ Res += (tmp * tmp.t()) / in.Ncols();
+ }
+ }
+ else
+ Res = cov(in.t());
+ return Res;
+ } //Matrix calc_corr
+
int do_work(int argc, char* argv[]) {
-
+
+ sttime = time(NULL);
Matrix MatrixData;
volume<float> Mean;
-
- {
+
+{
volume4D<float> RawData;
-
+ volume<float> theMask;
+ toc();
//read data
message("Reading data file " << (string)fnin.value() << " ... ");
read_volume4D(RawData,fnin.value());
message(" done" << endl);
-
+
Mean = meanvol(RawData);
-
+ toc();
+ message("Reading mask file " << (string)fnmask.value() << " ... ");
+ read_volume(theMask,fnmask.value());
+
memmsg(" Before reshape ");
- MatrixData = RawData.matrix();
+ MatrixData = RawData.matrix(theMask);
}
memmsg(" after reshape ");
message(" Data size " << MatrixData.Nrows() << " x " << MatrixData.Ncols() << endl);
-
+ toc();
memmsg(" before remmean_econ ");
remmean_econ(MatrixData);
memmsg("after remmean_econ / before remmean")
-
+ toc();
+
MatrixData = remmean(MatrixData);
-
memmsg(" after remmean ");
- message(" Mean size " << MatrixData.Nrows() << " x " << MatrixData.Ncols() << endl);
+ toc();
+ message(" Matrix size " << MatrixData.Nrows() << " x " << MatrixData.Ncols() << endl);
+
+ memmsg(" before calc_corr ");
+ toc();
+ calccorr(MatrixData,FALSE);
+ toc();
+ memmsg(" after calc_corr ");
+
+ memmsg(" before calc_corr 2");
+ toc();
+ calccorr(MatrixData,TRUE);
+ toc();
+ memmsg(" after calc_corr 2");
+
- message(" Mean size " << MatrixData.Nrows() << " x " << MatrixData.Ncols() << endl);
return 0;
}
@@ -117,6 +174,7 @@ int do_work(int argc, char* argv[]) {
// the help message is printed)
options.add(fnin);
+ options.add(fnmask);
options.add(help);
options.add(xdim);
options.add(ydim);