From 1edd63e1981ae94af675dc0d544c5062a46e8962 Mon Sep 17 00:00:00 2001
From: Pierre Kraemer <kraemer@modelis.u-strasbg.fr>
Date: Wed, 4 Jan 2012 17:26:15 +0100
Subject: [PATCH] particleCell2.5

---
 .../MovingObjects/particle_cell_2DandHalf.hpp | 40 ++++++++++---------
 include/Geometry/frame.h                      | 18 +++++----
 include/Geometry/transfo.h                    |  2 +-
 3 files changed, 33 insertions(+), 27 deletions(-)

diff --git a/include/Algo/MovingObjects/particle_cell_2DandHalf.hpp b/include/Algo/MovingObjects/particle_cell_2DandHalf.hpp
index 49d77cd3e..d9cf34fa9 100644
--- a/include/Algo/MovingObjects/particle_cell_2DandHalf.hpp
+++ b/include/Algo/MovingObjects/particle_cell_2DandHalf.hpp
@@ -1,5 +1,7 @@
 //#define DEBUG
 
+#include "Geometry/frame.h"
+
 namespace CGoGN
 {
 
@@ -42,7 +44,7 @@ typename PFP::VEC3 ParticleCell2DAndHalf<PFP>::pointInFace(Dart d)
 template <typename PFP>
 Geom::Orientation3D ParticleCell2DAndHalf<PFP>::getOrientationEdge(const VEC3& point, Dart d)
 {
-	const VEC3& endPoint = m_positions[m.phi2(d)];
+	const VEC3& endPoint = m_positions[m.phi1(d)];
 	const VEC3& vertexPoint = m_positions[d];
 
 	const VEC3& n1 = Algo::Geometry::faceNormal<PFP>(m, d, m_positions);
@@ -55,15 +57,15 @@ template <typename PFP>
 typename PFP::VEC3 ParticleCell2DAndHalf<PFP>::intersectLineEdge(const VEC3& pA, const VEC3& pB, Dart d)
 {
 	const VEC3& q1 = m_positions[d];
-	const VEC3& q2 = m_positions[m.phi2(d)];
+	const VEC3& q2 = m_positions[m.phi1(d)];
 	VEC3 Inter;
 
-	VEC3 n1 = Algo::Geometry::faceNormal<PFP>(m,d,m_positions);
+	VEC3 n1 = Algo::Geometry::faceNormal<PFP>(m, d, m_positions);
 	VEC3 n = (q2 - q1) ^ n1 ;
 
 	Geom::intersectionLinePlane(pA, pB - pA, q1, n, Inter) ;
 
-	Geom::Plane3D<float> pl = Algo::Geometry::facePlane<PFP>(m,d,m_positions);
+	Geom::Plane3D<float> pl = Algo::Geometry::facePlane<PFP>(m, d, m_positions);
 	pl.project(Inter);
 
 	return Inter;
@@ -172,32 +174,31 @@ void ParticleCell2DAndHalf<PFP>::edgeState(VEC3 current, Geom::Orientation3D sid
 	switch(sideOfEdge)
 	{
 		case Geom::UNDER :
+		{
 			d = m.phi1(d);
 			faceState(current);
 			return;
+		}
+
 		case Geom::OVER:
 		{
 			//transform the displacement into the new entered face
-			VEC3 displ = current-m_position;
-			VEC3 edge = Algo::Geometry::vectorOutOfDart<PFP>(m,m.phi2(d),m_positions);
-
-			VEC3 n1 = Algo::Geometry::faceNormal<PFP>(m,d,m_positions);
-			VEC3 n2 = Algo::Geometry::faceNormal<PFP>(m,m.phi2(d),m_positions);
-
-			float angle = acos(n1*n2);
+			VEC3 displ = current - m_position;
 
-			Geom::Matrix<4,4,float> mRot;
-			mRot.identity();
+			VEC3 n1 = Algo::Geometry::faceNormal<PFP>(m, d, m_positions);
+			VEC3 n2 = Algo::Geometry::faceNormal<PFP>(m, m.phi2(d), m_positions);
+			VEC3 axis = n1 ^ n2 ;
 
-			Geom::rotate(edge[0],edge[1],edge[2],angle,mRot);
+			float angle = Geom::angle(n1, n2) ;
 
-			displ = Geom::transform(displ,mRot);
-			current = m_position+displ;
+			displ = Geom::rotate(axis, angle, displ) ;
+			current = m_position + displ;
 
 			d = m.phi1(m.phi2(d));
 			faceState(current);
 			return;
 		}
+
 		default :
 			state = EDGE;
 	}
@@ -335,6 +336,7 @@ void ParticleCell2DAndHalf<PFP>::faceState(VEC3 current)
 	switch (getOrientationEdge(current, d))
 	{
 	case Geom::UNDER :
+		distance += (current - m_position).norm();
 		m_position = current;
 		state = FACE;
 		break;
@@ -347,12 +349,14 @@ void ParticleCell2DAndHalf<PFP>::faceState(VEC3 current)
 //
 //			vertexState(current);
 		}
-//		else
+		else
 		{
 // 			CGoGNout << "wsoe : " << wsoe << CGoGNendl;
 // 			CGoGNout << "current " << current << " " << m_position << CGoGNendl;
 // 			CGoGNout << "d " << d << "d1 " << m_positions[d] << " d2 " << m_positions[m.phi2(d)] << CGoGNendl;
-			m_position = intersectLineEdge(current, m_position, d);
+			VEC3 isect = intersectLineEdge(current, m_position, d);
+			distance += (isect - m_position).norm();
+			m_position = isect;
 // 			CGoGNout << " inter : " << m_position << CGoGNendl;
 			edgeState(current, Geom::OVER);
 		}
diff --git a/include/Geometry/frame.h b/include/Geometry/frame.h
index 7b5aafc6a..c8725c421 100644
--- a/include/Geometry/frame.h
+++ b/include/Geometry/frame.h
@@ -22,14 +22,16 @@
 *                                                                              *
 *******************************************************************************/
 
-#ifndef LOCALFRAME_H_
-#define LOCALFRAME_H_
+#ifndef _FRAME_H_
+#define _FRAME_H_
 
 #include <cmath>
 
-namespace CGoGN {
+namespace CGoGN
+{
 
-namespace Geom {
+namespace Geom
+{
 
 /**
  * Util for rotation of a 3D point (or vector) around a given line (going through the origin) and of a given angle
@@ -213,10 +215,10 @@ private : // private constants
 
 } ;
 
-} // Geom
+} // namespace Geom
 
-} // CGoGN
+} // namespace CGoGN
 
-#include "frame.hpp"
+#include "Geometry/frame.hpp"
 
-#endif /* LOCALFRAME_H_ */
+#endif /* _FRAME_H_ */
diff --git a/include/Geometry/transfo.h b/include/Geometry/transfo.h
index fed77e385..7879f9ef4 100644
--- a/include/Geometry/transfo.h
+++ b/include/Geometry/transfo.h
@@ -98,7 +98,7 @@ void rotate(Vector<3,T>& axis, T angle, Matrix<4,4,T>& mat);
  * @param mat the transformation matrix
  */
 template <typename T>
-Vector<3,T> transform(const Vector<3,T>& P,const Matrix<4,4,T>& mat);
+Vector<3,T> transform(const Vector<3,T>& P, const Matrix<4,4,T>& mat);
 
 } // namespace Geom
 
-- 
GitLab