diff --git a/vcg/complex/local_optimization/tri_edge_collapse_quadric.h b/vcg/complex/local_optimization/tri_edge_collapse_quadric.h
index b520c54f..5a4361e6 100644
--- a/vcg/complex/local_optimization/tri_edge_collapse_quadric.h
+++ b/vcg/complex/local_optimization/tri_edge_collapse_quadric.h
@@ -24,6 +24,9 @@
   History
 
 $Log: not supported by cvs2svn $
+Revision 1.13  2007/02/25 09:20:10  cignoni
+Added Rad to the NormalThr Option and removed a bug in multiple exectuion of non optimal simplification (missing an isD check)
+
 Revision 1.12  2007/01/19 09:13:14  cignoni
 Added Finalize() method to the interface, corrected minor bugs on border preserving and postsimplification cleanup
 Avoided double make_heap (it is done only in the local_optimization init)
@@ -168,7 +171,10 @@ public:
     typedef TriEdgeCollapseQuadricParameter QParameter;
     typedef HelperType QH;
 
-		static QParameter & Params(){static QParameter p; return p;}
+		static QParameter & Params(){
+			static QParameter p; 
+			return p;
+		}
 		enum Hint {
 			HNHasFFTopology       = 0x0001,  // La mesh arriva con la topologia ff gia'fatta
 			HNHasVFTopology       = 0x0002,  // La mesh arriva con la topologia bf gia'fatta
@@ -415,7 +421,8 @@ public:
 
     QuadricType qq=QH::Qd(v[0]);
     qq+=QH::Qd(v[1]);
-		double QuadErr = Params().ScaleFactor*qq.Apply(v[1]->P()); 
+		Point3d tpd=Point3d::Construct(v[1]->P());
+    double QuadErr = Params().ScaleFactor*qq.Apply(tpd); 
 
 		// All collapses involving triangles with quality larger than <QualityThr> has no penalty;
 		if(MinQual>Params().QualityThr) MinQual=Params().QualityThr;  
@@ -609,14 +616,17 @@ static void InitQuadric(TriMeshType &m)
 		q+=QH::Qd(v[1]);
 		
     Point3<QuadricType::ScalarType> x;
+		
     bool rt=q.Minimum(x);
 		if(!rt) { // if the computation of the minimum fails we choose between the two edge points and the middle one.
-			x.Import((v[0]->P()+v[1]->P())/2);
+			Point3<QuadricType::ScalarType> x0=Point3d::Construct(v[0]->P());
+		  Point3<QuadricType::ScalarType> x1=Point3d::Construct(v[1]->P());
+      x.Import((v[0]->P()+v[1]->P())/2);
 			double qvx=q.Apply(x);
-			double qv0=q.Apply(v[0]->P());
-			double qv1=q.Apply(v[1]->P());
-      if(qv0<qvx) x.Import(v[0]->P());
-			if(qv1<qvx && qv1<qv0) x.Import(v[1]->P());
+			double qv0=q.Apply(x0);
+			double qv1=q.Apply(x1);
+      if(qv0<qvx) x=x0;
+			if(qv1<qvx && qv1<qv0) x=x1;
 		}
 		
     return CoordType::Construct(x);
diff --git a/vcg/math/quadric.h b/vcg/math/quadric.h
index 0b740c3b..97d2ce6b 100644
--- a/vcg/math/quadric.h
+++ b/vcg/math/quadric.h
@@ -24,6 +24,9 @@
   History
 
 $Log: not supported by cvs2svn $
+Revision 1.7  2006/11/13 12:53:40  ponchio
+just added an #include <matrix33>
+
 Revision 1.6  2006/10/09 20:23:00  cignoni
 Added a minimum method that uses SVD. Unfortunately it is much much slower.
 
@@ -73,16 +76,16 @@ public:
 template< class PlaneType >
 	void ByPlane( const PlaneType & p )					// Init dato un piano
 	{
-		a[0] =  p.Direction()[0]*p.Direction()[0];	// a11
-		a[1] =  p.Direction()[1]*p.Direction()[0];	// a12 (=a21)
-		a[2] =  p.Direction()[2]*p.Direction()[0];	// a13 (=a31)
-		a[3] =  p.Direction()[1]*p.Direction()[1];	// a22
-		a[4] =  p.Direction()[2]*p.Direction()[1];	// a23 (=a32)
-		a[5] =  p.Direction()[2]*p.Direction()[2];	// a33
+		a[0] =  (ScalarType)p.Direction()[0]*p.Direction()[0];	// a11
+		a[1] =  (ScalarType)p.Direction()[1]*p.Direction()[0];	// a12 (=a21)
+		a[2] =  (ScalarType)p.Direction()[2]*p.Direction()[0];	// a13 (=a31)
+		a[3] =  (ScalarType)p.Direction()[1]*p.Direction()[1];	// a22
+		a[4] =  (ScalarType)p.Direction()[2]*p.Direction()[1];	// a23 (=a32)
+		a[5] =  (ScalarType)p.Direction()[2]*p.Direction()[2];	// a33
 		b[0] = (ScalarType)(-2.0)*p.Offset()*p.Direction()[0];
 		b[1] = (ScalarType)(-2.0)*p.Offset()*p.Direction()[1];
 		b[2] = (ScalarType)(-2.0)*p.Offset()*p.Direction()[2];
-		c    =  p.Offset()*p.Offset();
+		c    =  (ScalarType)p.Offset()*p.Offset();
 	}
 
 	void Zero()																// Azzera la quadrica
@@ -137,31 +140,32 @@ template <class ResultScalarType>
 	ResultScalarType Apply( const Point3<ResultScalarType> & p ) const	// Applica la quadrica al punto p
 	{
 		assert( IsValid() );
-
-	// Versione Lenta
 /*
+	// Versione Lenta
+
 		Point3d t;
 		t[0] = p[0]*a[0] + p[1]*a[1] + p[2]*a[2];
 		t[1] = p[0]*a[1] + p[1]*a[3] + p[2]*a[4];
 		t[2] = p[0]*a[2] + p[1]*a[4] + p[2]*a[5];
 		double k = b[0]*p[0] + b[1]*p[1] + b[2]*p[2];
-		double tp =t*p;
+		double tp = t*p ;
+		assert(tp+k+c >= 0);
 		return tp + k + c;
-	
-*/
-	/* Versione veloce */
+	*/
 
+	/* Versione veloce */
 		return ResultScalarType (
       p[0]*p[0]*a[0] + 2*p[0]*p[1]*a[1] + 2*p[0]*p[2]*a[2] + p[0]*b[0] 
 			               +   p[1]*p[1]*a[3] + 2*p[1]*p[2]*a[4] + p[1]*b[1]
 			                                  +   p[2]*p[2]*a[5] + p[2]*b[2]	+ c);
+
 	}
 
 // spostare..risolve un sistema 3x3
 template<class FLTYPE> 	
 bool Gauss33( FLTYPE x[], FLTYPE C[3][3+1] )
 {
-    const FLTYPE keps = (FLTYPE)1e-6;
+    const FLTYPE keps = (FLTYPE)1e-3;
     int i,j,k;
 
     FLTYPE eps;					// Determina valore cond.