diff --git a/pointmatcher/DataPointsFilters/SamplingSurfaceNormal.cpp b/pointmatcher/DataPointsFilters/SamplingSurfaceNormal.cpp
index 8f22a9b47..ee97de8c1 100644
--- a/pointmatcher/DataPointsFilters/SamplingSurfaceNormal.cpp
+++ b/pointmatcher/DataPointsFilters/SamplingSurfaceNormal.cpp
@@ -245,7 +245,7 @@ void SamplingSurfaceNormalDataPointsFilter<T>::fuseRange(
 	const Matrix NN = (d.colwise() - mean);
 
 	// compute covariance
-	const Matrix C(NN * NN.transpose());
+	const Matrix C((NN * NN.transpose()) / T(colCount));
 	Vector eigenVa = Vector::Identity(featDim-1, 1);
 	Matrix eigenVe = Matrix::Identity(featDim-1, featDim-1);
 	// Ensure that the matrix is suited for eigenvalues calculation
diff --git a/pointmatcher/DataPointsFilters/SurfaceNormal.cpp b/pointmatcher/DataPointsFilters/SurfaceNormal.cpp
index a836128e7..5ca4ea09a 100644
--- a/pointmatcher/DataPointsFilters/SurfaceNormal.cpp
+++ b/pointmatcher/DataPointsFilters/SurfaceNormal.cpp
@@ -184,7 +184,7 @@ void SurfaceNormalDataPointsFilter<T>::inPlaceFilter(
 		const Vector mean = d.rowwise().sum() / T(realKnn);
 		const Matrix NN = d.colwise() - mean;
 
-		const Matrix C(NN * NN.transpose());
+		const Matrix C((NN * NN.transpose()) / T(realKnn));
 		Vector eigenVa = Vector::Zero(featDim-1, 1);
 		Matrix eigenVe = Matrix::Zero(featDim-1, featDim-1);
 		// Ensure that the matrix is suited for eigenvalues calculation