Index: cpp/frams/model/similarity/SVD/matrix_tools.cpp
===================================================================
--- cpp/frams/model/similarity/SVD/matrix_tools.cpp	(revision 816)
+++ cpp/frams/model/similarity/SVD/matrix_tools.cpp	(revision 817)
@@ -84,5 +84,5 @@
 }
 
-double *Transpose(double *&A, int arow, int acol)
+double *Transpose(double *&A, int arow, int acol, double *&toDel, int delSize)
 {
 	double *result = Create(acol * arow);
@@ -94,9 +94,36 @@
 		}
 
+	if (delSize != 0)
+		delete[] toDel;
+	
 	return result;
-
-}
-
-/** Computes the SVD of the nSize x nSize distance matrix
+}
+
+//Weighted centring of a matrix.
+//https://github.com/vegandevs/vegan/blob/master/src/goffactor.c
+void wcentre(double *x, double *w, int *nr, int *nc)
+{
+	int i, j, ij;
+	double sw, swx;
+
+	for (i = 0, sw = 0.0; i < (*nr); i++)
+		sw += w[i];
+
+	for (j = 0; j < (*nc) ; j++)
+	{
+		for (i = 0, swx = 0.0, ij = (*nr)*j; i < (*nr); i++, ij++)
+		{
+			swx += w[i] * x[ij];
+		}
+		swx /= sw;
+		for (i = 0,  ij = (*nr)*j; i < (*nr); i++, ij++)
+		{
+			x[ij] -= swx;
+			x[ij] *= sqrt(w[i]);
+		}
+	}
+}
+
+/** Computes the weighted MDS of the nSize x nSize distance matrix
 		@param vdEigenvalues [OUT] Vector of doubles. On return holds the eigenvalues of the
 		decomposed distance matrix (or rather, to be strict, of the matrix of scalar products
@@ -106,86 +133,47 @@
 		@param pDistances [IN] matrix of distances between parts.
 		@param Coordinates [OUT] array of three dimensional coordinates obtained from SVD of pDistances matrix.
+  		@param weights [IN] vector of row weights.
 		*/
-void MatrixTools::SVD(std::vector<double> &vdEigenvalues, int nSize, double *pDistances, Pt3D *&Coordinates)
-{
-	// compute squares of elements of this array
-	// compute the matrix B that is the parameter of SVD
-	double *B = Create(nSize * nSize);
-	{
-		// use additional scope to delete temporary matrices
-		double *Ones, *Eye, *Z, *D;
-
-		D = Create(nSize * nSize);
-		D = Power(pDistances, nSize, nSize, 2.0, D, nSize);
-
-		Ones = Create(nSize * nSize);
-		for (int i = 0; i < nSize; i++)
-			for (int j = 0; j < nSize; j++)
-			{
-			Ones[i * nSize + j] = 1;
-			}
-
-		Eye = Create(nSize * nSize);
-		for (int i = 0; i < nSize; i++)
-		{
-			for (int j = 0; j < nSize; j++)
-			{
-				if (i == j)
-				{
-					Eye[i * nSize + j] = 1;
-				}
-				else
-				{
-					Eye[i * nSize + j] = 0;
-				}
-			}
-		}
-
-		Z = Create(nSize * nSize);
-		for (int i = 0; i < nSize; i++)
-		{
-			for (int j = 0; j < nSize; j++)
-			{
-				Z[i * nSize + j] = 1.0 / ((double)nSize) * Ones[i * nSize + j];
-			}
-		}
-
-		for (int i = 0; i < nSize; i++)
-		{
-			for (int j = 0; j < nSize; j++)
-			{
-				Z[i * nSize + j] = Eye[i * nSize + j] - Z[i * nSize + j];
-			}
-		}
-
-		for (int i = 0; i < nSize; i++)
-		{
-			for (int j = 0; j < nSize; j++)
-			{
-				B[i * nSize + j] = Z[i * nSize + j] * -0.5;
-			}
-		}
-
-		B = Multiply(B, D, nSize, nSize, nSize, B, nSize);
-		B = Multiply(B, Z, nSize, nSize, nSize, B, nSize);
-
-		delete[] Ones;
-		delete[] Eye;
-		delete[] Z;
-		delete[] D;
-	}
+void MatrixTools::weightedMDS(std::vector<double> &vdEigenvalues, int nSize, double *pDistances, Pt3D *&Coordinates, double *weights)
+{
+	// compute the matrix D that is the parameter of SVD
+	double *D = Create(nSize * nSize);
+	D = Power(pDistances, nSize, nSize, 2.0, D, nSize);
+
+	for (int i = 0; i < 2; i++)
+	{
+		wcentre(D, weights, &nSize, &nSize);
+		D = Transpose(D, nSize, nSize, D, nSize);
+	}
+		
+	for (int i = 0; i < nSize; i++)
+		for (int j = 0; j < nSize; j++)
+		{
+			D[i * nSize + j] *= -0.5;
+		}
 
 	double *Eigenvalues = Create(nSize);
 	double *S = Create(nSize * nSize);
 
-	// call SVD function
+	// call the SVD function
 	double *Vt = Create(nSize * nSize);
 	size_t astep = nSize * sizeof(double);
-	Lapack::JacobiSVD(B, astep, Eigenvalues, Vt, astep, nSize, nSize, nSize);
-
-	double *W = Transpose(Vt, nSize, nSize);
-
-	delete[] B;
+	Lapack::JacobiSVD(D, astep, Eigenvalues, Vt, astep, nSize, nSize, nSize);
+
+	double *W = Transpose(Vt, nSize, nSize, W, 0);
+
+	delete[] D;
 	delete[] Vt;
+	
+	// deweight
+	double row_weight = 1;
+	for (int i = 0; i < nSize; i++)
+	{
+		row_weight = weights[i];
+		for (int j = 0; j < nSize; j++)
+		{
+			W[i * nSize + j] /= sqrt(row_weight);
+		}
+	}
 
 	for (int i = 0; i < nSize; i++)
Index: cpp/frams/model/similarity/SVD/matrix_tools.h
===================================================================
--- cpp/frams/model/similarity/SVD/matrix_tools.h	(revision 816)
+++ cpp/frams/model/similarity/SVD/matrix_tools.h	(revision 817)
@@ -12,5 +12,5 @@
 {
 public:
-	static void SVD(std::vector<double> &vdEigenvalues, int nSize, double *pDistances, Pt3D *&Coordinates);
+	static void weightedMDS(std::vector<double> &vdEigenvalues, int nSize, double *pDistances, Pt3D *&Coordinates, double *weights);
 };
 
Index: cpp/frams/model/similarity/simil_model.cpp
===================================================================
--- cpp/frams/model/similarity/simil_model.cpp	(revision 816)
+++ cpp/frams/model/similarity/simil_model.cpp	(revision 817)
@@ -35,5 +35,5 @@
 
 static ParamEntry MSparam_tab[] = {
-		{ "Creature: Similarity", 1, 6, "ModelSimilarity", "Evaluates morphological dissimilarity. More information:\nhttp://www.framsticks.com/bib/Komosinski-et-al-2001\nhttp://www.framsticks.com/bib/Komosinski-and-Kubiak-2011\nhttp://www.framsticks.com/bib/Komosinski-2016", },
+		{ "Creature: Similarity", 1, 7, "ModelSimilarity", "Evaluates morphological dissimilarity. More information:\nhttp://www.framsticks.com/bib/Komosinski-et-al-2001\nhttp://www.framsticks.com/bib/Komosinski-and-Kubiak-2011\nhttp://www.framsticks.com/bib/Komosinski-2016", },
 		{ "simil_parts", 0, 0, "Weight of parts count", "f 0 100 0", FIELD(m_adFactors[0]), "Differing number of parts is also handled by the 'part degree' similarity component.", },
 		{ "simil_partdeg", 0, 0, "Weight of parts' degree", "f 0 100 1", FIELD(m_adFactors[1]), "", },
@@ -41,4 +41,5 @@
 		{ "simil_partgeom", 0, 0, "Weight of parts' geometric distances", "f 0 100 0", FIELD(m_adFactors[3]), "", },
 		{ "simil_fixedZaxis", 0, 0, "Fix 'z' (vertical) axis?", "d 0 1 0", FIELD(fixedZaxis), "", },
+		{ "simil_weightedMDS", 0, 0, "Should weighted MDS be used?", "d 0 1 0", FIELD(wMDS), "", },
 		{ "evaluateDistance", 0, PARAM_DONTSAVE | PARAM_USERHIDDEN, "evaluate model dissimilarity", "p f(oGeno,oGeno)", PROCEDURE(p_evaldistance), "Calculates dissimilarity between two models created from Geno objects.", },
 		{ 0, },
@@ -75,4 +76,7 @@
 	isFuzzy = 0;
 	fuzzyDepth = 10;
+	
+	//Determines whether weighted MDS should be used.
+	wMDS = 0;
 }
 
@@ -1981,5 +1985,5 @@
 		int nSize = m_Mod[iMod]->getPartCount();
 
-		double *pDistances = (double *)malloc(nSize * nSize * sizeof(double));
+		double *pDistances =  new double[nSize * nSize];
 
 		for (int i = 0; i < nSize; i++)
@@ -1993,4 +1997,21 @@
 		Pt3D P1Pos, P2Pos; // positions of parts 
 		double dDistance; // the distance between Parts
+		
+		double *weights = new double[nSize];
+		for (int i = 0; i < nSize; i++)
+		{
+			if (wMDS==1)
+				weights[i] = 0;
+			else
+				weights[i] = 1;
+		}
+				
+		if (wMDS==1)
+			for (int i = 0; i < pModel->getJointCount(); i++)
+			{
+				weights[pModel->getJoint(i)->p1_refno]++;
+				weights[pModel->getJoint(i)->p2_refno]++; 
+			}
+		
 		for (iP1 = 0; iP1 < pModel->getPartCount(); iP1++)
 		{
@@ -2020,5 +2041,5 @@
 		} // for (iP1)
 
-		MatrixTools::SVD(vEigenvalues, nSize, pDistances, m_aPositions[iMod]);
+		MatrixTools::weightedMDS(vEigenvalues, nSize, pDistances, m_aPositions[iMod], weights);
 		if (fixedZaxis == 1) //restore the original vertical coordinate of each Part
 		{
@@ -2029,5 +2050,6 @@
 		}
 
-		free(pDistances);
+		delete[] pDistances;
+		delete[] weights;
 	}
 
Index: cpp/frams/model/similarity/simil_model.h
===================================================================
--- cpp/frams/model/similarity/simil_model.h	(revision 816)
+++ cpp/frams/model/similarity/simil_model.h	(revision 817)
@@ -88,4 +88,7 @@
 	int fuzzyDepth;
 	int isFuzzy;
+        
+	//For wMDS = 1 weighted MDS with vertex degrees as weights is used for the alignment.
+	int wMDS;
 
 	/// Interface to local parameters
