Skip to content
Snippets Groups Projects
Commit 913c5bdf authored by Christian Beckmann's avatar Christian Beckmann
Browse files

update instacorr

parent 8de55d55
No related branches found
No related tags found
No related merge requests found
...@@ -120,7 +120,7 @@ namespace Melodic{ ...@@ -120,7 +120,7 @@ namespace Melodic{
tmpTC = tmpData * insta_maps.t(); tmpTC = tmpData * insta_maps.t();
if(opts.insta_fn.value().length()>0){ if(opts.insta_fn.value().length()>0){
dbgmsg(string("BEGIN: INSTACORR") << endl);
if(opts.insta_idx.value()<1 || opts.insta_idx.value()>tmpTC.Ncols()){ if(opts.insta_idx.value()<1 || opts.insta_idx.value()>tmpTC.Ncols()){
cerr << "ERROR:: INSTACORR index is wrong \n\n"; cerr << "ERROR:: INSTACORR index is wrong \n\n";
exit(2); exit(2);
...@@ -129,25 +129,31 @@ namespace Melodic{ ...@@ -129,25 +129,31 @@ namespace Melodic{
Matrix tmpRef = tmpTC.Column(opts.insta_idx.value()); Matrix tmpRef = tmpTC.Column(opts.insta_idx.value());
if(opts.insta_idx.value()>1){ if(opts.insta_idx.value()>1){
// swap columns // swap columns
dbgmsg(string("INSTACORR: swap columns") << endl);
tmpTC.Column(opts.insta_idx.value()) << tmpTC.Column(1); tmpTC.Column(opts.insta_idx.value()) << tmpTC.Column(1);
tmpTC.Column(1) << tmpRef; tmpTC.Column(1) << tmpRef;
} }
if(opts.insta_partial.value() && tmpTC.Ncols()>1){ if(opts.insta_partial.value() && tmpTC.Ncols()>1){
// partal correlations // partal correlations
dbgmsg(string("INSTACORR: partial analysis") << endl);
Matrix tmpConf = tmpTC.Columns(2,tmpTC.Ncols()); Matrix tmpConf = tmpTC.Columns(2,tmpTC.Ncols());
tmpData -= tmpConf * (pinv(tmpConf) * tmpData); tmpData -= tmpConf * (pinv(tmpConf) * tmpData);
tmpRef -= tmpConf * (pinv(tmpConf) * tmpRef); tmpRef -= tmpConf * (pinv(tmpConf) * tmpRef);
} }
if(opts.insta_varnorm.value()){ if(opts.insta_varnorm.value()){
dbgmsg(string("INSTACORR: varnorm") << endl);
Matrix vscales = pow(stdev(tmpData,1),-1); Matrix vscales = pow(stdev(tmpData,1),-1);
varnorm(tmpData,vscales); varnorm(tmpData,vscales);
varnorm(tmpRef,pow(stdev(tmpRef,1),-1)); varnorm(tmpRef,pow(stdev(tmpRef,1),-1));
} }
// Shur product // Shur product
dbgmsg(string("INSTACORR: SP") << endl);
SP4(tmpData,tmpRef); SP4(tmpData,tmpRef);
dbgmsg(string("END: INSTACORR") << endl);
} }
//END INSTACORRS //END INSTACORRS
......
...@@ -66,7 +66,7 @@ ...@@ -66,7 +66,7 @@
namespace Melodic{ namespace Melodic{
const string version = "3.14"; const string version = "3.141";
// The two strings below specify the title and example usage that is // The two strings below specify the title and example usage that is
// printed out as the help or usage message // printed out as the help or usage message
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment