diff --git a/meldata.cc b/meldata.cc
index 86dba71471d80ec2032241ca2548e53be35a13a6..65c6dbddca6a0f59abc34ebe8047644d0abc8630 100644
--- a/meldata.cc
+++ b/meldata.cc
@@ -168,21 +168,21 @@ namespace Melodic{
       add_Smodes(tmpS2);
     }
 
-		//add GLM OLS fit
-		if(Tdes.Storage()){
-			Matrix alltcs = Tmodes.at(0).Column(1);
-			for(int ctr=1; ctr < (int)Tmodes.size();ctr++)
-				alltcs|=Tmodes.at(ctr).Column(1);
-			if((alltcs.Nrows()==Tdes.Nrows())&&(Tdes.Nrows()>Tdes.Ncols()))
-				glmT.olsfit(alltcs,Tdes,Tcon);
-		}
-		if(Sdes.Storage()){
-			Matrix alltcs = Smodes.at(0);
-			for(int ctr=1; ctr < (int)Smodes.size();ctr++)
-				alltcs|=Smodes.at(ctr);
-			if((alltcs.Nrows()==Sdes.Nrows())&&(Sdes.Nrows()>Sdes.Ncols()&&alltcs.Nrows()>2))
-				glmS.olsfit(alltcs,Sdes,Scon);
-		}
+	//add GLM OLS fit
+	if(Tdes.Storage()){
+		Matrix alltcs = Tmodes.at(0).Column(1);
+		for(int ctr=1; ctr < (int)Tmodes.size();ctr++)
+			alltcs|=Tmodes.at(ctr).Column(1);
+		if((alltcs.Nrows()==Tdes.Nrows())&&(Tdes.Nrows()>Tdes.Ncols()))
+			glmT.olsfit(alltcs,Tdes,Tcon);
+	}
+	if(Sdes.Storage()){
+		Matrix alltcs = Smodes.at(0);
+		for(int ctr=1; ctr < (int)Smodes.size();ctr++)
+			alltcs|=Smodes.at(ctr);
+		if((alltcs.Nrows()==Sdes.Nrows())&&(Sdes.Nrows()>Sdes.Ncols()&&alltcs.Nrows()>2))
+			glmS.olsfit(alltcs,Sdes,Scon);
+	}
   }
 
   void MelodicData::setup()
@@ -256,9 +256,17 @@ namespace Melodic{
 			PPCA=tmpPPCA;
   		}
     	order = opts.pca_dim.value();
-    	calc_white(pcaE, pcaD, order, whiteMatrix, dewhiteMatrix);
 
-    	if(opts.debug.value()){
+		if (opts.paradigmfname.value().length()>0){
+			Matrix tmpPscales;
+			tmpPscales = param.t() * alldat;
+			paramS = stdev(tmpPscales.t());
+		    
+    		calc_white(pcaE, pcaD, order, param, paramS, whiteMatrix, dewhiteMatrix);
+		}else
+    		calc_white(pcaE, pcaD, order, whiteMatrix, dewhiteMatrix);
+
+		if(opts.debug.value()){
 			outMsize("pcaE",pcaE); saveascii(pcaE,"pcaE");
 			outMsize("pcaD",pcaD); saveascii(pcaD,"pcaD");
 			outMsize("AdjEV",AdjEV); saveascii(AdjEV,"AdjEV");
@@ -332,16 +340,16 @@ namespace Melodic{
     //create mask
     create_mask(Mask);
 
-		//setup background image
-		if(opts.bgimage.value()>""){
-			read_volume(background,opts.bgimage.value());
-      if(!samesize(Mean,background)){
-        cerr << "ERROR:: background image and data have different dimensions  \n\n";
-        exit(2);
-      }
-		}else{
-			background = Mean;
-		}
+	//setup background image
+	if(opts.bgimage.value()>""){
+		read_volume(background,opts.bgimage.value());
+    	if(!samesize(Mean,background)){
+        	cerr << "ERROR:: background image and data have different dimensions  \n\n";
+        	exit(2);
+    	}
+	}else{
+		background = Mean;
+	}
 		
     if(!samesize(Mean,Mask)){
       cerr << "ERROR:: mask and data have different dimensions  \n\n";
@@ -360,28 +368,46 @@ namespace Melodic{
     double tmptime = time(NULL);
     srand((unsigned int) tmptime);
 
-		//read in post-proc design matrices etc
-		if(opts.fn_Tdesign.value().length()>0)
-			Tdes = read_ascii_matrix(opts.fn_Tdesign.value());
-		if(opts.fn_Sdesign.value().length()>0)
-			Sdes = read_ascii_matrix(opts.fn_Sdesign.value());
-		if(opts.fn_Tcon.value().length()>0)
-			Tcon = read_ascii_matrix(opts.fn_Tcon.value());
-		if(opts.fn_Scon.value().length()>0)
-			Scon = read_ascii_matrix(opts.fn_Scon.value());
-		if(opts.fn_TconF.value().length()>0)
-			TconF = read_ascii_matrix(opts.fn_TconF.value());
-		if(opts.fn_SconF.value().length()>0)
-			SconF = read_ascii_matrix(opts.fn_SconF.value());
+	if(opts.paradigmfname.value().length()>0){
+		message("  Use columns in " << opts.paradigmfname.value() 
+	      << " for PCA initialisation" <<endl);
+		param = read_ascii_matrix(opts.paradigmfname.value());
 		
-		if(numfiles>1 && Sdes.Storage() == 0){
-			Sdes = ones(numfiles,1);
-			if(Scon.Storage() == 0){
-				Scon = ones(1,1);
-				Scon &= -1*Scon;
-			}
+	    Matrix tmpPU, tmpPV;
+		DiagonalMatrix tmpPD;
+		SVD(param, tmpPD, tmpPU, tmpPV);
+		param = tmpPU;
+		
+		opts.pca_dim.set_T(std::max(opts.pca_dim.value(), param.Ncols()+3));		
+		if(opts.debug.value()){
+			outMsize("Paradigm",param); saveascii(param,"param");
+			outMsize("ParadigmS",paramS); saveascii(param,"paramS");
 		}
-		Tdes = remmean(Tdes,1);
+		opts.guessfname.set_T(opts.paradigmfname.value());
+	}
+
+	//read in post-proc design matrices etc
+	if(opts.fn_Tdesign.value().length()>0)
+		Tdes = read_ascii_matrix(opts.fn_Tdesign.value());
+	if(opts.fn_Sdesign.value().length()>0)
+		Sdes = read_ascii_matrix(opts.fn_Sdesign.value());
+	if(opts.fn_Tcon.value().length()>0)
+		Tcon = read_ascii_matrix(opts.fn_Tcon.value());
+	if(opts.fn_Scon.value().length()>0)
+		Scon = read_ascii_matrix(opts.fn_Scon.value());
+	if(opts.fn_TconF.value().length()>0)
+		TconF = read_ascii_matrix(opts.fn_TconF.value());
+	if(opts.fn_SconF.value().length()>0)
+		SconF = read_ascii_matrix(opts.fn_SconF.value());
+		
+	if(numfiles>1 && Sdes.Storage() == 0){
+ 		Sdes = ones(numfiles,1);
+		if(Scon.Storage() == 0){
+			Scon = ones(1,1);
+			Scon &= -1*Scon;
+		}
+	}
+	Tdes = remmean(Tdes,1);
   }
 
   void MelodicData::save()
@@ -585,7 +611,8 @@ namespace Melodic{
     RXweight = tmpRX.matrix(Mask);
   } 
 
-  void MelodicData::est_smoothness(){
+  void MelodicData::est_smoothness()
+  {
     if(Resels == 0){
       string SM_path = opts.binpath + "smoothest";
       string Mask_fname = logger.appendDir("mask");
@@ -620,7 +647,8 @@ namespace Melodic{
     }
   }
 
-  unsigned long MelodicData::standardise(volume<float>& mask, volume4D<float>& R){
+  unsigned long MelodicData::standardise(volume<float>& mask, volume4D<float>& R)
+  {
     	unsigned long count = 0;
     	int M=R.tsize();
     
@@ -809,7 +837,8 @@ namespace Melodic{
     }
   } //void create_mask()
 
-  void MelodicData::sort(){
+  void MelodicData::sort()
+  {
     int numComp = mixMatrix.Ncols(), numVox = IC.Ncols(), 
         numTime = mixMatrix.Nrows(), i,j;
 
@@ -891,11 +920,13 @@ namespace Melodic{
     unmixMatrix = pinv(mixMatrix);
   }
 
-  void MelodicData::reregress(){
+  void MelodicData::reregress()
+  {
 	if((numfiles > 1)){
 	   
 	}
   }
+  
   void MelodicData::status(const string &txt)
   {
     cout << "MelodicData Object " << txt << endl;
diff --git a/meldata.h b/meldata.h
index 37eb6d5cb9e66cd8ec353a60def9ca46a3ffe159..c4b911bb163721df0887812546b925df26a87908 100644
--- a/meldata.h
+++ b/meldata.h
@@ -94,6 +94,12 @@ namespace Melodic{
 
       void set_TSmode();
 
+      inline Matrix& get_param() {return param;} 
+      inline void set_param(Matrix& Arg) {param = Arg;}	
+
+      inline Matrix& get_paramS() {return paramS;} 
+      inline void set_paramS(Matrix& Arg) {paramS = Arg;}	
+
       inline Matrix& get_white() {return whiteMatrix;}
       inline void set_white(Matrix& Arg) {whiteMatrix = Arg;}
       
@@ -194,7 +200,7 @@ namespace Melodic{
       volumeinfo tempInfo;
       vector<Matrix> DWM, WM;
 			basicGLM glmT, glmS;
-			Matrix Tdes, Tcon, TconF, Sdes, Scon, SconF;	
+			Matrix Tdes, Tcon, TconF, Sdes, Scon, SconF, param, paramS;	
 			RowVector explained_var;		
 
     private:
diff --git a/melhlprfns.cc b/melhlprfns.cc
index 1a81a796435b80c1f74454315a941c4116768a42..d160abcd32f0486eac96a8560f86e3bd48003f1d 100644
--- a/melhlprfns.cc
+++ b/melhlprfns.cc
@@ -169,7 +169,7 @@ namespace Melodic{
     return Res;
   }  //Matrix SP
 
-	Matrix corrcoef(const Matrix& in1, const Matrix& in2){
+  Matrix corrcoef(const Matrix& in1, const Matrix& in2){
 		Matrix tmp = in1;
 		tmp |= in2;
 		Matrix out;
@@ -219,6 +219,50 @@ 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)
+  {
+
+//	tmpD2= tmpD | tmpPD.AsRow().Columns(tmpPE.Ncols()-param.Ncols()+1,tmpPE.Ncols());
+//  cerr << tmpPD.AsRow().Columns(tmpPE.Ncols()-param.Ncols()+1,tmpPE.Ncols()) << endl;
+
+//    
+
+//	Matrix tmpPE;
+//	tmpPE = SP(param,ones(param.Nrows(),1)*pow(stdev(param,1)*std::sqrt((float)param.Nrows()),-1));
+
+//	RE |= tmpPE;
+//	RowVector tmpD2;
+//	tmpD2 = tmpD | stdev(param,1).AsRow()*std::sqrt((float)param.Nrows());
+//	RD << abs(diag(tmpD2.t()));
+//    RD << RD.SymSubMatrix(N-dim+1,N);
+
+	Matrix RE;
+    DiagonalMatrix RD;
+    int N = tmpE.Ncols();
+    dim = std::min(dim,N);
+
+//	cerr << stdev(param) << endl;
+    RE = tmpE.Columns(std::min(N-dim+1+param.Ncols(),N-2),N);
+	RE |= param;
+
+//	cerr << paramS.Nrows() << " x " << paramS.Ncols() << endl;
+//	cerr << paramS << endl;
+	RowVector tmpD2;
+	tmpD2 = tmpD | pow(paramS,2).AsRow();
+    RD << abs(diag(tmpD2.t()));
+
+//	cerr << " " <<tmpD2.Ncols() << " " << N << " " << dim << endl;
+    RD << RD.SymSubMatrix(N-dim+1+param.Ncols(),N+param.Ncols());    
+
+    float res = 1.0;    
+    white = sqrt(abs(RD.i()))*RE.t();
+    dewhite = RE*sqrt(RD);
+
+    if(dim < N)
+      res = PercEV(dim);
+    return res;
+  }  //Matrix calc_white
+
   float calc_white(const Matrix& tmpE, const RowVector& tmpD, const RowVector& PercEV, int dim, Matrix& white, Matrix& dewhite)
   {
     Matrix RE;
@@ -245,6 +289,13 @@ namespace Melodic{
     tmp2 = calc_white(tmpE,tmpD, tmp, dim, white, dewhite); 
   }  //Matrix calc_white
 
+  void calc_white(const Matrix& tmpE, const RowVector& tmpD, int dim, Matrix& param, Matrix& paramS, Matrix& white, Matrix& dewhite)
+  {
+    RowVector tmp(tmpE.Ncols());
+    float tmp2;
+    tmp2 = calc_white(tmpE,tmpD, tmp, dim, param, paramS, white, dewhite); 
+  }  //Matrix calc_white
+
   void calc_white(const Matrix& Corr, int dim, Matrix& white, Matrix& dewhite)
   {
     Matrix RE;
@@ -257,6 +308,7 @@ namespace Melodic{
     calc_white(RE,tmp2, dim, white, dewhite); 
   }  //Matrix calc_white
   
+ 
   void std_pca(const Matrix& Mat, const Matrix& weights, Matrix& Corr, Matrix& evecs, RowVector& evals)
   {
     if(weights.Storage()>0)
diff --git a/melhlprfns.h b/melhlprfns.h
index 1d51cab11464d2d45d08b977969b9ccd759324e4..791b4d27c1f1a1a8a39a09bd07c56b61f8b7eab3 100644
--- a/melhlprfns.h
+++ b/melhlprfns.h
@@ -29,7 +29,7 @@ namespace Melodic{
   Matrix convert_to_pbsc(Matrix& Mat);
 
   RowVector varnorm(Matrix& in, int dim = 30, float level = 1.6);
-	void varnorm(Matrix& in, const RowVector& vars);
+       void varnorm(Matrix& in, const RowVector& vars);
   RowVector varnorm(Matrix& in, Matrix& Corr, int dim = 30, float level = 1.6);
 
   Matrix SP2(const Matrix& in, const Matrix& weights, bool econ = 0);
@@ -37,11 +37,13 @@ namespace Melodic{
   RowVector Feta(int n1,int n2);
   RowVector cumsum(const RowVector& Inp);
 
-	Matrix corrcoef(const Matrix& in1, const Matrix& in2);
+  Matrix corrcoef(const Matrix& in1, const Matrix& in2);
   Matrix calc_corr(const Matrix& in, bool econ = 0);
   Matrix calc_corr(const Matrix& in, const Matrix& weights, bool econ = 0);
 
+  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);
+  void calc_white(const Matrix& tmpE, const RowVector& tmpD, int dim, Matrix& param, Matrix& paramS, Matrix& white, Matrix& dewhite);
   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);
   
diff --git a/melodic.cc b/melodic.cc
index 4935a2fa53ec6d188bbde0cbda2d0ed95660df86..e46501e43c6db0470b2082a76f8a5157b930f24b 100644
--- a/melodic.cc
+++ b/melodic.cc
@@ -96,8 +96,8 @@ int main(int argc, char *argv[]){
 				no_conv = icaobj.no_convergence;
 
 				opts.maxNumItt.set_T(500);
-				if((opts.approach.value()=="symm")&&
-	   		(retry > std::min(opts.retrystep,3))){
+	  if((opts.approach.value()=="symm")&&(retry > std::min(opts.retrystep,3)))
+	  {
 	  			if(no_conv){
 	    			retry++;
 	    			opts.approach.set_T("defl");
@@ -166,7 +166,7 @@ int main(int argc, char *argv[]){
 
 				message("finished!" << endl << endl);
       } 
-			else { 
+	  else { 
 				message(endl <<"No convergence -- giving up " << endl << endl);
       }	     
     }
@@ -277,8 +277,7 @@ void mmonly(Log& logger, MelodicOptions& opts,
     mmres = mmall(logger,opts,melodat,report,pmaps);
 	}
 
-Matrix mmall(Log& logger, MelodicOptions& opts,
-	MelodicData& melodat, MelodicReport& report, Matrix& pmaps){
+Matrix mmall(Log& logger, MelodicOptions& opts,MelodicData& melodat, MelodicReport& report, Matrix& pmaps){
   
   Matrix mmpars(5*melodat.get_IC().Nrows(),5);
   mmpars = 0;
@@ -307,7 +306,7 @@ Matrix mmall(Log& logger, MelodicOptions& opts,
       ICmap = SP(melodat.get_IC().Row(ctr),diagvals(ctr)*melodat.get_stdNoisei());
     }else{
     	ICmap = melodat.get_IC().Row(ctr);
-		}
+	}
     string wherelog;
     if(opts.genreport.value())
       wherelog = report.getDir();
diff --git a/meloptions.h b/meloptions.h
index 08ff4213fb759debd19a52eca923e7ea4fe49d27..ae6bffcf5bccdf30d0eb653a045d336140bccec9 100644
--- a/meloptions.h
+++ b/meloptions.h
@@ -155,11 +155,11 @@ class MelodicOptions {
 
  inline MelodicOptions::MelodicOptions() :
    logdir(string("-o,--outdir"), string("log.txt"),
-	  string("output directory name\n"), 
-	  false, requires_argument),
+	   string("output directory name\n"), 
+	   false, requires_argument),
    inputfname(string("-i,--in"), std::vector<string>(),
-	      string("input file names (either single file name or comma-separated list or text file)"), 
-	      true, requires_argument),
+	   string("input file names (either single file name or comma-separated list or text file)"), 
+	   true, requires_argument),
    outputfname(string("-O,--out"), string("melodic"),
 	   string("output file name"), 
 	   false, requires_argument,false),
@@ -187,12 +187,12 @@ class MelodicOptions {
    joined_whiten(string("--sep_whiten"), true,
 	   string("switch on separate whitening"), 
 	   false, no_argument),
-	 joined_vn(string("--sep_vn"), true,
-		 string("switch off joined variance nomalisation"), 
-		 false, no_argument),
-	 vn_level(string("--vn_level"), float(2.3),
-		 string("variance nomalisation threshold level (Z> value is ignored)"), 
-		 false, requires_argument, false),
+   joined_vn(string("--sep_vn"), true,
+   	   string("switch off joined variance nomalisation"), 
+       false, no_argument),
+   vn_level(string("--vn_level"), float(2.3),
+	   string("variance nomalisation threshold level (Z> value is ignored)"), 
+	   false, requires_argument, false),
    numICs(string("-n,--numICs"), -1,
 	   string("numer of IC's to extract (for deflation approach)"), 
 	   false, requires_argument),
@@ -220,9 +220,9 @@ class MelodicOptions {
    epsilon(string("--eps,--epsilon"), 0.0005,
 	   string("minimum error change"), 
 	   false, requires_argument),
- 	 epsilonS(string("--epsS,--epsilonS"), 0.03,
-		   string("minimum error change for rank-1 approximation in TICA"), 
-		   false, requires_argument),
+   epsilonS(string("--epsS,--epsilonS"), 0.03,
+	   string("minimum error change for rank-1 approximation in TICA"), 
+	   false, requires_argument),
    maxNumItt(string("--maxit"),  500,
 	   string("\tmaximum number of iterations before restart"), 
 	   false, requires_argument),
@@ -247,54 +247,54 @@ class MelodicOptions {
    smodename(string("--smode"),  string(""),
 	   string("\tmatrix of session modes for report generation"), 
 	   false, requires_argument),
-	 filter(string("-f,--filter"),  string(""),
-		 string("component numbers to remove\n "), 
-		 false, requires_argument),
+   filter(string("-f,--filter"),  string(""),
+	   string("component numbers to remove\n "), 
+	   false, requires_argument),
    genreport(string("--report"), false,
 	   string("generate Melodic web report"), 
 	   false, no_argument),
-	 guireport(string("--guireport"), string(""),
-		   string("modify report for GUI use"), 
-		   false, requires_argument, false),
-	 bgimage(string("--bgimage"),  string(""),
-		 string("specify background image for report (default: mean image)\n "), 
-		 false, requires_argument),
+   guireport(string("--guireport"), string(""),
+	   string("modify report for GUI use"), 
+	   false, requires_argument, false),
+   bgimage(string("--bgimage"),  string(""),
+	   string("specify background image for report (default: mean image)\n "), 
+	   false, requires_argument),
    tr(string("--tr"),  0.0,
 	   string("\tTR in seconds"), 
 	   false, requires_argument),
-	 logPower(string("--logPower"),  false,
-		 string("calculate log of power for frequency spectrum\n"), 
-		 false, no_argument),
-	 addsigchng(string("--sigchng"),  false,
-		 string("add signal change estimates to report pages\n"), 
-		 false, no_argument, false),
-	 allPPCA(string("--allPPCA"),  false,
-			 string("add all PPCA plots\n"), 
-			 false, no_argument, false),
-	 varplots(string("--varplots"),  false,
-		 string("add std error envelopes to time course plots\n"), 
-		 false, no_argument, false),
-	 varvals(string("--varvals"),  false,
-		 string("add rank1 values after plots\n"), 
-		 false, no_argument, false),
-	 fn_Tdesign(string("--Tdes"), string(""),
+   logPower(string("--logPower"),  false,
+	   string("calculate log of power for frequency spectrum\n"), 
+	   false, no_argument),
+   addsigchng(string("--sigchng"),  false,
+	   string("add signal change estimates to report pages\n"), 
+       false, no_argument, false),
+   allPPCA(string("--allPPCA"),  false,
+	   string("add all PPCA plots\n"), 
+	   false, no_argument, false),
+   varplots(string("--varplots"),  false,
+	   string("add std error envelopes to time course plots\n"), 
+	   false, no_argument, false),
+   varvals(string("--varvals"),  false,
+	   string("add rank1 values after plots\n"), 
+	   false, no_argument, false),
+   fn_Tdesign(string("--Tdes"), string(""),
 	   string("        design matrix across time-domain"),
-		 false, requires_argument),
-	 fn_Tcon(string("--Tcon"), string(""),
-		 string("        t-contrast matrix across time-domain"),
-		 false, requires_argument),
+	   false, requires_argument),
+   fn_Tcon(string("--Tcon"), string(""),
+       string("        t-contrast matrix across time-domain"),
+	   false, requires_argument),
    fn_TconF(string("--Tconf"), string(""),
-		 string("        F-contrast matrix across time-domain"),
-		 false, requires_argument, false),
+	   string("        F-contrast matrix across time-domain"),
+	   false, requires_argument, false),
    fn_Sdesign(string("--Sdes"), string(""),
-		 string("        design matrix across subject-domain"),
-		 false, requires_argument),
-	 fn_Scon(string("--Scon"), string(""),
-		 string("        t-contrast matrix across subject-domain"),
-		 false, requires_argument),	
-	 fn_SconF(string("--Sconf"), string(""),
-		 string("        F-contrast matrix across subject-domain"),
-		 false, requires_argument,false),	
+	   string("        design matrix across subject-domain"),
+	   false, requires_argument),
+   fn_Scon(string("--Scon"), string(""),
+	   string("        t-contrast matrix across subject-domain"),
+	   false, requires_argument),	
+   fn_SconF(string("--Sconf"), string(""),
+	   string("        F-contrast matrix across subject-domain"),
+	   false, requires_argument,false),	
    output_all(string("--Oall"),  false,
 	   string("        output everything"), 
 	   false, no_argument),
@@ -335,7 +335,7 @@ class MelodicOptions {
 	   string("file name of guess of mixing matrix"), 
 	   false, requires_argument, false),
    paradigmfname(string("-p,--paradigm"),  string(""),
-	   string("file name of FEAT paradigm file"), 
+	   string("file name of FEAT paradigm file (design.mat)"), 
 	   false, requires_argument, false),
    dummy(string("--dummy"),  0,
 	   string("number of dummy volumes"), 
@@ -355,21 +355,21 @@ class MelodicOptions {
    remove_meanvol(string("--keep_meanvol"), true,
 	   string("do not subtract mean volume"), 
 	   false, no_argument, false),
-   remove_meantc(string("--keep_meantc"), true,
-		   string("do not remove mean time course"), 
-		   false, no_argument, false),
+   remove_meantc(string("--remove_meantc"), false,
+	   string("remove mean time course"), 
+	   false, no_argument, false),
    remove_endslices(string("--remEndslices"),  false,
 	   string("delete end slices (motion correction artefacts)"), 
 	   false, no_argument,false),
    rescale_nht(string("--rescale_nht"),  true,
 	   string("switch off map rescaling after mixture-modelling"), 
 	   false, no_argument,false),
-   guess_remderiv(string("--remderiv"),  false,
+   guess_remderiv(string("--remove_deriv"),  false,
 	   string("removes every second entry in paradigm file (EV derivatives)"), 
-	   false, no_argument, false),
+	   false, no_argument),
    temporal(string("--temporal"),  false,
 	   string("perform temporal ICA"), 
-	    false, no_argument, false),
+	   false, no_argument, false),
    retrystep(3),
    options(title, usageexmpl)
    {
@@ -409,19 +409,19 @@ class MelodicOptions {
 	    options.add(filter);
 	    options.add(genreport);
 	    options.add(guireport);
-			options.add(bgimage);
+		options.add(bgimage);
 	    options.add(tr);
 	    options.add(logPower);
 	    options.add(addsigchng);
 	    options.add(allPPCA);
 	    options.add(varplots);
 	    options.add(varvals);
-			options.add(fn_Tdesign);
-			options.add(fn_Tcon);
-			options.add(fn_TconF);
-			options.add(fn_Sdesign);
-			options.add(fn_Scon);
-			options.add(fn_SconF);
+		options.add(fn_Tdesign);
+		options.add(fn_Tcon);
+		options.add(fn_TconF);
+		options.add(fn_Sdesign);
+		options.add(fn_Scon);
+		options.add(fn_SconF);
 	    options.add(output_all);
 	    options.add(output_unmix);
 	    options.add(output_MMstats);
diff --git a/melpca.cc b/melpca.cc
index 4c1edb9206a989768ffc0fa71c49d709f8a459c8..3c703ee7aada9e918f1a32e891beff32b2895dc3 100644
--- a/melpca.cc
+++ b/melpca.cc
@@ -32,7 +32,17 @@ namespace Melodic{
     Matrix tmpE;
     RowVector tmpD, AdjEV, PercEV;
    
-    std_pca(in,weights,Corr,tmpE,tmpD); 
+	if(opts.paradigmfname.value().length()>0)
+	{
+		basicGLM tmpglm;
+		tmpglm.olsfit(in,melodat.get_param(),Identity(melodat.get_param().Ncols()));
+		std_pca(tmpglm.get_residu(),weights,Corr,tmpE,tmpD);
+//		std_pca(in,weights,Corr,tmpE,tmpD,melodat.get_param());
+	}
+	else{
+ 		std_pca(in,weights,Corr,tmpE,tmpD); 
+	}	
+	
     if(opts.tsmooth.value()){
       message(endl << "  temporal smoothing of Eigenvectors " << endl);
       tmpE=smoothColumns(tmpE);
@@ -80,8 +90,13 @@ namespace Melodic{
     Matrix tmpDeWhite;
 
     float varp = 1.0;
-    varp = calc_white(melodat.get_pcaE(),melodat.get_pcaD(), 
-    	melodat.get_EVP() ,opts.pca_dim.value(),tmpWhite,tmpDeWhite);
+
+	if (opts.paradigmfname.value().length()>0)
+	    varp = calc_white(melodat.get_pcaE(),melodat.get_pcaD(), 
+    		melodat.get_EVP(),opts.pca_dim.value(),melodat.get_param(),melodat.get_paramS(),tmpWhite,tmpDeWhite);
+	else
+		varp = calc_white(melodat.get_pcaE(),melodat.get_pcaD(), 
+    		melodat.get_EVP(),opts.pca_dim.value(),tmpWhite,tmpDeWhite);
 
     melodat.set_white(tmpWhite);
     melodat.set_dewhite(tmpDeWhite);