Création d'un compte pour un collaborateur extérieur au laboratoire depuis l'intranet ICube : https://intranet.icube.unistra.fr/fr/labs/member/profile

surfaceDeformation.cpp 19.3 KB
Newer Older
1
2
#include "surfaceDeformation.h"

3
4
5
#include "Algo/Selection/raySelector.h"
#include "Algo/Selection/collector.h"

6
7
8
#include "Algo/Geometry/normal.h"
#include "Algo/Geometry/laplacian.h"

9
10
11
#include <QKeyEvent>
#include <QMouseEvent>

12
13
namespace CGoGN
{
14

15
namespace SCHNApps
16
{
17

18
19
20
21
22
23
PerMapParameterSet::PerMapParameterSet(MapHandlerGen* m) :
	mh(m),
	verticesSelectionMode(LOCKED),
	nlContext(NULL)
{
	QString positionName;
24

25
	QString vec3TypeName = QString::fromStdString(nameOfType(PFP2::VEC3()));
26

27
28
	const AttributeHash& attribs = mh->getAttributesList(VERTEX);
	for(AttributeHash::const_iterator i = attribs.constBegin(); i != attribs.constEnd(); ++i)
29
	{
30
		if(i.value() == vec3TypeName)
31
		{
32
33
			if(positionName != "position")	// try to select an attribute named "position"
				positionName = i.key();		// or anything else if not found
34
35
		}
	}
36
	positionAttribute = mh->getAttribute<PFP2::VEC3, VERTEX>(positionName);
37

38
	PFP2::MAP* map = static_cast<MapHandler<PFP2>*>(mh)->getMap();
39
40
41
	lockingMarker = new CellMarker<VERTEX>(*map);
	handleMarker = new CellMarker<VERTEX>(*map);

42
	positionInit = mh->getAttribute<PFP2::VEC3, VERTEX>("positionInit", false) ;
43
	if(!positionInit.isValid())
44
		positionInit = mh->addAttribute<PFP2::VEC3, VERTEX>("positionInit", false) ;
45

46
47
48
	vIndex = mh->getAttribute<unsigned int, VERTEX>("vIndex", false);
	if(!vIndex.isValid())
		vIndex = mh->addAttribute<unsigned int, VERTEX>("vIndex", false);
49

50
51
52
//	edgeWeight = mh->getAttribute<PFP2::REAL, EDGE>("edgeWeight", false);
//	if(!edgeWeight.isValid())
//		edgeWeight = mh->addAttribute<PFP2::REAL, EDGE>("edgeWeight", false);
53

54
55
56
//	vertexArea = mh->getAttribute<PFP2::REAL, VERTEX>("vertexArea", false);
//	if(!vertexArea.isValid())
//		vertexArea = mh->addAttribute<PFP2::REAL, VERTEX>("vertexArea", false);
57

58
	diffCoord = mh->getAttribute<PFP2::VEC3, VERTEX>("diffCoord", false);
59
	if(!diffCoord.isValid())
60
		diffCoord = mh->addAttribute<PFP2::VEC3, VERTEX>("diffCoord", false);
61

62
	vertexRotationMatrix = mh->getAttribute<Eigen_Matrix3f, VERTEX>("vertexRotationMatrix", false) ;
63
	if(!vertexRotationMatrix.isValid())
64
		vertexRotationMatrix = mh->addAttribute<Eigen_Matrix3f, VERTEX>("vertexRotationMatrix", false);
65

66
	rotatedDiffCoord = mh->getAttribute<PFP2::VEC3, VERTEX>("rotatedDiffCoord", false) ;
67
	if(!rotatedDiffCoord.isValid())
68
		rotatedDiffCoord = mh->addAttribute<PFP2::VEC3, VERTEX>("rotatedDiffCoord", false);
69

70
71
	initParameters();
}
72

73
74
75
76
77
78
79
PerMapParameterSet::~PerMapParameterSet()
{
	delete lockingMarker;
	delete handleMarker;
	nlDeleteContext(nlContext);
}

80
81
82
83
84
void PerMapParameterSet::initParameters()
{
	if(positionAttribute.isValid())
	{
		PFP2::MAP* map = static_cast<MapHandler<PFP2>*>(mh)->getMap();
85

86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
		map->copyAttribute(positionInit, positionAttribute) ;

		nb_vertices = static_cast<AttribMap*>(map)->computeIndexCells<VERTEX>(vIndex);

//		Algo::Surface::Geometry::computeCotanWeightEdges<PFP2>(*map, positionAttribute, edgeWeight) ;
//		Algo::Surface::Geometry::computeVoronoiAreaVertices<PFP2>(*map, positionAttribute, vertexArea) ;
//		Algo::Surface::Geometry::computeLaplacianCotanVertices<PFP2>(*map, edgeWeight, vertexArea, positionAttribute, diffCoord) ;
		Algo::Surface::Geometry::computeLaplacianTopoVertices<PFP2>(*map, positionAttribute, diffCoord) ;

		for(unsigned int i = vertexRotationMatrix.begin(); i != vertexRotationMatrix.end() ; vertexRotationMatrix.next(i))
			vertexRotationMatrix[i] = Eigen::Matrix3f::Identity();

		if(nlContext)
			nlDeleteContext(nlContext);
		nlContext = nlNewContext();
		nlSolverParameteri(NL_NB_VARIABLES, nb_vertices);
		nlSolverParameteri(NL_LEAST_SQUARES, NL_TRUE);
		nlSolverParameteri(NL_SOLVER, NL_CHOLMOD_EXT);
	}
105
106
107
}


108
109
bool SurfaceDeformationPlugin::enable()
{
110
	m_dockTab = new SurfaceDeformationDockTab(m_window, this);
111
112
	addTabInDock(m_dockTab, "Surface Deformation");

113
114
	m_drawer = new Utils::Drawer();
	registerShader(m_drawer->getShader());
115
116
117
118
119

	connect(m_window, SIGNAL(viewAndPluginLinked(View*, Plugin*)), this, SLOT(viewLinked(View*, Plugin*)));
	connect(m_window, SIGNAL(viewAndPluginUnlinked(View*, Plugin*)), this, SLOT(viewUnlinked(View*, Plugin*)));
	connect(m_window, SIGNAL(currentViewChanged(View*)), this, SLOT(currentViewChanged(View*)));

120
121
122
	return true;
}

123
124
void SurfaceDeformationPlugin::disable()
{
125
	delete m_drawer;
126
127
128
129
130
131
}

void SurfaceDeformationPlugin::redraw(View* view)
{
	if(selecting)
	{
132
		glDisable(GL_LIGHTING) ;
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
		m_drawer->newList(GL_COMPILE_AND_EXECUTE) ;
		m_drawer->lineWidth(2.0f) ;
		m_drawer->begin(GL_LINES) ;
		m_drawer->color3f(0.0f, 0.0f, 1.0f) ;
		m_drawer->vertex(selectionCenter) ;
		m_drawer->vertex(selectionCenter + selectionRadius * PFP2::VEC3(1,0,0)) ;
		m_drawer->vertex(selectionCenter) ;
		m_drawer->vertex(selectionCenter + selectionRadius * PFP2::VEC3(-1,0,0)) ;
		m_drawer->vertex(selectionCenter) ;
		m_drawer->vertex(selectionCenter + selectionRadius * PFP2::VEC3(0,1,0)) ;
		m_drawer->vertex(selectionCenter) ;
		m_drawer->vertex(selectionCenter + selectionRadius * PFP2::VEC3(0,-1,0)) ;
		m_drawer->vertex(selectionCenter) ;
		m_drawer->vertex(selectionCenter + selectionRadius * PFP2::VEC3(0,0,1)) ;
		m_drawer->vertex(selectionCenter) ;
		m_drawer->vertex(selectionCenter + selectionRadius * PFP2::VEC3(0,0,-1)) ;
		m_drawer->end() ;
		m_drawer->endList() ;
151
	}
152
153
154
155
156
157
158
159
160

	ParameterSet* params = h_viewParams[view];
	MapHandlerGen* mh = params->selectedMap;
	if(mh)
	{
		PerMapParameterSet* perMap = params->perMap[mh->getName()];

		if(!perMap->locked_vertices.empty() || !perMap->handle_vertices.empty())
		{
161
			glDisable(GL_LIGHTING) ;
162
163
164
165
166
167
168
169
170
171
172
173
174
175
			m_drawer->newList(GL_COMPILE_AND_EXECUTE) ;
			m_drawer->pointSize(4.0f) ;
			m_drawer->begin(GL_POINTS) ;
			m_drawer->color3f(1.0f, 0.0f, 0.0f) ;
			for(unsigned int i = 0; i < perMap->locked_vertices.size(); ++i)
			{
				if (!perMap->handleMarker->isMarked(perMap->locked_vertices[i]))
					m_drawer->vertex(perMap->positionAttribute[perMap->locked_vertices[i]]) ;
			}
			m_drawer->color3f(0.0f, 0.0f, 1.0f) ;
			for(unsigned int i = 0; i < perMap->handle_vertices.size(); ++i)
				m_drawer->vertex(perMap->positionAttribute[perMap->handle_vertices[i]]) ;
			m_drawer->end() ;
			m_drawer->endList() ;
176
		}
177
178
179
	}
}

180
181
182
183
184
void SurfaceDeformationPlugin::keyPress(View* view, QKeyEvent* event)
{
	if(event->key() == Qt::Key_Shift)
	{
		view->setMouseTracking(true);
185
		selecting = true;
186
187
		view->updateGL();
	}
188
	else if(event->key() == Qt::Key_R)
189
190
191
	{
		ParameterSet* params = h_viewParams[view];
		MapHandlerGen* map = params->selectedMap;
192
193
194
		if(map)
		{
			asRigidAsPossible(view, map);
195
196
			PerMapParameterSet* perMap = params->perMap[map->getName()];
			params->selectedMap->notifyAttributeModification(perMap->positionAttribute);
197
198
			view->updateGL();
		}
199
	}
200
201
202
203
204
205
206
}

void SurfaceDeformationPlugin::keyRelease(View* view, QKeyEvent* event)
{
	if(event->key() == Qt::Key_Shift)
	{
		view->setMouseTracking(false);
207
		selecting = false;
208
209
210
211
212
213
		view->updateGL();
	}
}

void SurfaceDeformationPlugin::mousePress(View* view, QMouseEvent* event)
{
214
	if(event->button() == Qt::LeftButton && selecting)
215
216
	{
		ParameterSet* params = h_viewParams[view];
217
		PerMapParameterSet* perMap = params->getCurrentMapParameterSet();
218

219
220
221
222
223
224
		if(perMap)
		{
			QPoint pixel(event->x(), event->y());
			qglviewer::Vec orig;
			qglviewer::Vec dir;
			view->camera()->convertClickToLine(pixel, orig, dir);
225

226
227
			PFP2::VEC3 rayA(orig.x, orig.y, orig.z);
			PFP2::VEC3 AB(dir.x, dir.y, dir.z);
228

229
230
231
			Dart d ;
			PFP2::MAP* map = static_cast<MapHandler<PFP2>*>(params->selectedMap)->getMap();
			Algo::Selection::vertexRaySelection<PFP2>(*map, perMap->positionAttribute, rayA, AB, d) ;
232

233
			if(d != NIL)
234
			{
235
236
237
238
239
				Algo::Surface::Selection::Collector_WithinSphere<PFP2> neigh(*map, perMap->positionAttribute, selectionRadius) ;
				neigh.collectAll(d) ;
				const std::vector<Dart>& insideV = neigh.getInsideVertices() ;

				for(unsigned int i = 0; i < insideV.size(); ++i)
240
				{
241
242
243
244
245
246
247
					unsigned int v = map->getEmbedding<VERTEX>(insideV[i]) ;
					if (!perMap->lockingMarker->isMarked(v))
					{
						perMap->locked_vertices.push_back(v) ;
						perMap->lockingMarker->mark(v);
					}
					if(perMap->verticesSelectionMode == HANDLE)
248
					{
249
250
251
252
253
						if(!perMap->handleMarker->isMarked(v))
						{
							perMap->handle_vertices.push_back(v) ;
							perMap->handleMarker->mark(v);
						}
254
255
					}
				}
256

257
258
				nlMakeCurrent(perMap->nlContext) ;
				nlReset(NL_FALSE) ;
259

260
261
				view->updateGL() ;
			}
262
263
		}
	}
264
265
266
267
268
269
270
	else if(event->button() == Qt::RightButton && event->modifiers() & Qt::ShiftModifier)
	{
		view->setMouseTracking(false) ;

		ParameterSet* params = h_viewParams[view];
		PerMapParameterSet* perMap = params->getCurrentMapParameterSet();

271
		if(perMap)
272
		{
273
274
			selecting = false ;
			dragging = true ;
275

276
277
278
279
280
281
282
283
284
285
286
287
			dragZ = 0;
			for(unsigned int i = 0; i < perMap->handle_vertices.size(); ++i)
			{
				const PFP2::VEC3& p = perMap->positionAttribute[perMap->handle_vertices[i]] ;
				qglviewer::Vec q = view->camera()->projectedCoordinatesOf(qglviewer::Vec(p[0],p[1],p[2]));
				dragZ += q.z ;
			}
			dragZ /= perMap->handle_vertices.size() ;

			qglviewer::Vec p(event->x(), event->y(), dragZ);
			dragPrevious = view->camera()->unprojectedCoordinatesOf(p);
		}
288
	}
289
290
291
292
293
294
}

void SurfaceDeformationPlugin::mouseRelease(View* view, QMouseEvent* event)
{
	if(event->button() == Qt::RightButton)
	{
295
		dragging = false ;
296
297
298
299
300
301
	}
}

void SurfaceDeformationPlugin::mouseMove(View* view, QMouseEvent* event)
{
	ParameterSet* params = h_viewParams[view];
302
	PerMapParameterSet* perMap = params->getCurrentMapParameterSet();
303

304
	if(perMap)
305
	{
306
307
308
309
		if(dragging)
		{
			qglviewer::Vec p(event->x(), event->y(), dragZ);
			qglviewer::Vec q = view->camera()->unprojectedCoordinatesOf(p);
310

311
312
313
314
			qglviewer::Vec vec = q - dragPrevious;
			PFP2::VEC3 t(vec.x, vec.y, vec.z);
			for(unsigned int i = 0; i < perMap->handle_vertices.size(); ++i)
				perMap->positionAttribute[perMap->handle_vertices[i]] += t ;
315

316
			dragPrevious = q ;
317

318
319
//			matchDiffCoord(view, map);
			asRigidAsPossible(view, params->selectedMap);
320

321
			params->selectedMap->notifyAttributeModification(perMap->positionAttribute);
322

323
324
325
			view->updateGL();
		}
		else if(selecting)
326
		{
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
			QPoint pixel(event->x(), event->y());
			qglviewer::Vec orig;
			qglviewer::Vec dir;
			view->camera()->convertClickToLine(pixel, orig, dir);

			PFP2::VEC3 rayA(orig.x, orig.y, orig.z);
			PFP2::VEC3 AB(dir.x, dir.y, dir.z);

			Dart d ;
			PFP2::MAP* map = static_cast<MapHandler<PFP2>*>(params->selectedMap)->getMap();
			Algo::Selection::vertexRaySelection<PFP2>(*map, perMap->positionAttribute, rayA, AB, d) ;
			if(d != NIL)
			{
				selectionCenter = perMap->positionAttribute[d] ;
				view->updateGL() ;
			}
343
344
345
346
		}
	}
}

347
348
349
350
351
352
353
354
355
356
357
358
void SurfaceDeformationPlugin::wheelEvent(View* view, QWheelEvent* event)
{
	if(selecting)
	{
		if(event->delta() > 0)
			selectionRadius *= 0.9f ;
		else
			selectionRadius *= 1.1f ;
		view->updateGL() ;
	}
}

359
360
361
362
363
364
365
366
367
void SurfaceDeformationPlugin::viewLinked(View* view, Plugin* plugin)
{
	if(plugin == this)
	{
		ParameterSet* params = new ParameterSet();
		h_viewParams.insert(view, params);
		const QList<MapHandlerGen*>& maps = view->getLinkedMaps();
		foreach(MapHandlerGen* map, maps)
		{
368
			PerMapParameterSet* p = new PerMapParameterSet(map);
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
			params->perMap.insert(map->getName(), p);
		}
		if (!maps.empty())
			changeSelectedMap(view, maps[0]);

		connect(view, SIGNAL(mapLinked(MapHandlerGen*)), this, SLOT(mapLinked(MapHandlerGen*)));
		connect(view, SIGNAL(mapUnlinked(MapHandlerGen*)), this, SLOT(mapUnlinked(MapHandlerGen*)));

		if(view->isCurrentView())
			m_dockTab->refreshUI(params);
	}
}

void SurfaceDeformationPlugin::viewUnlinked(View* view, Plugin* plugin)
{
	if(plugin == this)
	{
386
387
388
389
390
391
392
		ParameterSet* params = h_viewParams[view];
		QHash<QString, PerMapParameterSet*>::const_iterator i = params->perMap.constBegin();
		while (i != params->perMap.constEnd())
		{
			delete i.value();
			++i;
		}
393
394
395
396
397
398
399
400
401
402
		h_viewParams.remove(view);

		disconnect(view, SIGNAL(mapLinked(MapHandlerGen*)), this, SLOT(mapLinked(MapHandlerGen*)));
		disconnect(view, SIGNAL(mapUnlinked(MapHandlerGen*)), this, SLOT(mapUnlinked(MapHandlerGen*)));
	}
}

void SurfaceDeformationPlugin::currentViewChanged(View* view)
{
	if(isLinkedToView(view))
403
404
405
406
407
	{
		ParameterSet* params = h_viewParams[view];
		changeSelectedMap(view, params->selectedMap);
		m_dockTab->refreshUI(params);
	}
408
409
410
411
412
413
414
415
}

void SurfaceDeformationPlugin::mapLinked(MapHandlerGen* m)
{
	View* view = static_cast<View*>(QObject::sender());
	assert(isLinkedToView(view));

	ParameterSet* params = h_viewParams[view];
416
	PerMapParameterSet* p = new PerMapParameterSet(m);
417
418
419
420
421
422
423
424
425
426
427
428
429
	params->perMap.insert(m->getName(), p);
	if(params->selectedMap == NULL || params->perMap.count() == 1)
		changeSelectedMap(view, m);
	else
		m_dockTab->refreshUI(params);
}

void SurfaceDeformationPlugin::mapUnlinked(MapHandlerGen* m)
{
	View* view = static_cast<View*>(QObject::sender());
	assert(isLinkedToView(view));

	ParameterSet* params = h_viewParams[view];
430
	delete params->perMap[m->getName()];
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
	params->perMap.remove(m->getName());

	if(params->selectedMap == m)
	{
		if(!params->perMap.empty())
			changeSelectedMap(view, m_window->getMap(params->perMap.begin().key()));
		else
			changeSelectedMap(view, NULL);
	}
	else
		m_dockTab->refreshUI(params);
}

void SurfaceDeformationPlugin::changeSelectedMap(View* view, MapHandlerGen* map)
{
	ParameterSet* params = h_viewParams[view];

	MapHandlerGen* prev = params->selectedMap;
	params->selectedMap = map;

	if(view->isCurrentView())
	{
		if(prev)
454
			disconnect(prev, SIGNAL(attributeAdded(unsigned int, const QString&)), m_dockTab, SLOT(addAttributeToList(unsigned int, const QString&)));
455
		if(map)
456
		{
457
			connect(map, SIGNAL(attributeAdded(unsigned int, const QString&)), m_dockTab, SLOT(addAttributeToList(unsigned int, const QString&)));
458
459
			selectionRadius = map->getBBdiagSize() / 50.0;
		}
460

461
462
463
464
465
466
467
		m_dockTab->refreshUI(params);
	}
}

void SurfaceDeformationPlugin::changePositionAttribute(View* view, MapHandlerGen* map, VertexAttribute<PFP2::VEC3> attribute)
{
	ParameterSet* params = h_viewParams[view];
468
469
470
	PerMapParameterSet* perMap = params->perMap[map->getName()];
	perMap->positionAttribute = attribute;
	perMap->initParameters();
471
472
473

	if(view->isCurrentView())
		m_dockTab->refreshUI(params);
474
475
476
477
478
479
480
481
482
}

void SurfaceDeformationPlugin::changeVerticesSelectionMode(View* view, MapHandlerGen* map, SelectionMode m)
{
	ParameterSet* params = h_viewParams[view];
	params->perMap[map->getName()]->verticesSelectionMode = m;

	if(view->isCurrentView())
		m_dockTab->refreshUI(params);
483
484
}

485
486
487
488
489
void SurfaceDeformationPlugin::matchDiffCoord(View* view, MapHandlerGen* mh)
{
	PFP2::MAP* map = static_cast<MapHandler<PFP2>*>(mh)->getMap();
	PerMapParameterSet* perMap = h_viewParams[view]->perMap[mh->getName()];

490
491
492
	nlMakeCurrent(perMap->nlContext);
	if(nlGetCurrentState() == NL_STATE_INITIAL)
		nlBegin(NL_SYSTEM) ;
493
494
	for(int coord = 0; coord < 3; ++coord)
	{
495
496
		LinearSolving::setupVariables<PFP2>(*map, perMap->vIndex, *perMap->lockingMarker, perMap->positionAttribute, coord) ;
		nlBegin(NL_MATRIX) ;
497
498
		LinearSolving::addRowsRHS_Laplacian_Topo<PFP2>(*map, perMap->vIndex, perMap->diffCoord, coord) ;
//		LinearSolving::addRowsRHS_Laplacian_Cotan<PFP2>(*map, perMap->vIndex, perMap->edgeWeight, perMap->vertexArea, perMap->diffCoord, coord) ;
499
500
501
502
503
		nlEnd(NL_MATRIX) ;
		nlEnd(NL_SYSTEM) ;
		nlSolve() ;
		LinearSolving::getResult<PFP2>(*map, perMap->vIndex, perMap->positionAttribute, coord) ;
		nlReset(NL_TRUE) ;
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
	}
}

void SurfaceDeformationPlugin::asRigidAsPossible(View* view, MapHandlerGen* mh)
{
	PFP2::MAP* map = static_cast<MapHandler<PFP2>*>(mh)->getMap();
	PerMapParameterSet* perMap = h_viewParams[view]->perMap[mh->getName()];

	CellMarkerNoUnmark<VERTEX> m(*map) ;

	for(Dart d = map->begin(); d != map->end(); map->next(d))
	{
		if(!m.isMarked(d))
		{
			m.mark(d) ;

			Eigen::Matrix3f cov = Eigen::Matrix3f::Zero() ;
			PFP2::VEC3 p = perMap->positionAttribute[d] ;
			PFP2::VEC3 pInit = perMap->positionInit[d] ;
523
//			PFP2::REAL area = perMap->vertexArea[d] ;
524
525
526
527
528
529
530
531
			Dart it = d ;
			do
			{
				Dart neigh = map->phi1(it) ;
				PFP2::VEC3 v = perMap->positionAttribute[neigh] - p ;
				PFP2::VEC3 vv = perMap->positionInit[neigh] - pInit ;
				for(unsigned int i = 0; i < 3; ++i)
					for(unsigned int j = 0; j < 3; ++j)
532
						cov(i,j) += v[i] * vv[j];// * perMap->edgeWeight[it] / area ;
533
534
535
536
537
538
539
				Dart dboundary = map->phi_1(it) ;
				if(map->phi2(dboundary) == dboundary)
				{
					v = perMap->positionAttribute[dboundary] - p ;
					vv = perMap->positionInit[dboundary] - p ;
					for(unsigned int i = 0; i < 3; ++i)
						for(unsigned int j = 0; j < 3; ++j)
540
							cov(i,j) += v[i] * vv[j];// * perMap->edgeWeight[dboundary] / area ;
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
				}
				it = map->alpha1(it) ;
			} while(it != d) ;

			Eigen::JacobiSVD<Eigen::Matrix3f> svd(cov, Eigen::ComputeFullU | Eigen::ComputeFullV) ;
			Eigen::Matrix3f R = svd.matrixU() * svd.matrixV().transpose() ;

			if(R.determinant() < 0)
			{
				Eigen::Matrix3f U = svd.matrixU() ;
				for(unsigned int i = 0; i < 3; ++i)
					U(i,2) *= -1 ;
				R = U * svd.matrixV().transpose() ;
			}

			perMap->vertexRotationMatrix[d] = R ;
		}
	}

	for(Dart d = map->begin(); d != map->end(); map->next(d))
	{
		if(m.isMarked(d))
		{
			m.unmark(d) ;

			unsigned int degree = 0 ;
			Eigen::Matrix3f r = Eigen::Matrix3f::Zero() ;
			Dart it = d ;
			do
			{
				r += perMap->vertexRotationMatrix[map->phi1(it)] ;
				++degree ;
				Dart dboundary = map->phi_1(it) ;
				if(map->phi2(dboundary) == dboundary)
				{
					r += perMap->vertexRotationMatrix[dboundary] ;
					++degree ;
				}
				it = map->alpha1(it) ;
			} while(it != d) ;
			r += perMap->vertexRotationMatrix[d] ;
			r /= degree + 1 ;
			PFP2::VEC3& dc = perMap->diffCoord[d] ;
			Eigen::Vector3f rdc(dc[0], dc[1], dc[2]) ;
			rdc = r * rdc ;
			perMap->rotatedDiffCoord[d] = PFP2::VEC3(rdc[0], rdc[1], rdc[2]) ;

//			Eigen::Vector3f rdc = Eigen::Vector3f::Zero() ;
//			Dart it = d ;
//			PFP::REAL vArea = perMap->vertexArea[d] ;
//			Eigen::Matrix3f& rotM = perMap->vertexRotationMatrix[d] ;
//			PFP::REAL val = 0 ;
//			do
//			{
//				Dart ddn = map->phi1(it) ;
//				PFP::REAL w = perMap->edgeWeight[it] / vArea ;
//				PFP::VEC3 vv = (perMap->positionInit[ddn] - perMap->positionInit[it]) * w ;
//				val += w ;
//				Eigen::Matrix3f r = 0.5f * (perMap->vertexRotationMatrix[ddn] + rotM) ;
//				Eigen::Vector3f vvr = r * Eigen::Vector3f(vv[0], vv[1], vv[2]) ;
//				rdc += vvr ;
//				Dart dboundary = map->phi_1(it) ;
//				if(map->phi2(dboundary) == dboundary)
//				{
//					w = perMap->edgeWeight[dboundary] / vArea ;
//					vv = (perMap->positionInit[dboundary] - perMap->positionInit[it]) * w ;
//					val += w ;
//					r = 0.5f * (perMap->vertexRotationMatrix[dboundary] + rotM) ;
//					vvr = r * Eigen::Vector3f(vv[0], vv[1], vv[2]) ;
//					rdc += vvr ;
//				}
//				it = map->alpha1(it) ;
//			} while(it != d) ;
//			rdc /= val ;
//			perMap->rotatedDiffCoord[d] = PFP::VEC3(rdc[0], rdc[1], rdc[2]) ;
		}
	}

619
620
621
	nlMakeCurrent(perMap->nlContext);
	if(nlGetCurrentState() == NL_STATE_INITIAL)
		nlBegin(NL_SYSTEM) ;
622
623
	for(int coord = 0; coord < 3; ++coord)
	{
624
625
		LinearSolving::setupVariables<PFP2>(*map, perMap->vIndex, *perMap->lockingMarker, perMap->positionAttribute, coord) ;
		nlBegin(NL_MATRIX) ;
626
627
//		LinearSolving::addRowsRHS_Laplacian_Cotan<PFP2>(*map, perMap->vIndex, perMap->edgeWeight, perMap->vertexArea, perMap->rotatedDiffCoord, coord) ;
		LinearSolving::addRowsRHS_Laplacian_Topo<PFP2>(*map, perMap->vIndex, perMap->rotatedDiffCoord, coord) ;
628
629
630
631
632
		nlEnd(NL_MATRIX) ;
		nlEnd(NL_SYSTEM) ;
		nlSolve() ;
		LinearSolving::getResult<PFP2>(*map, perMap->vIndex, perMap->positionAttribute, coord) ;
		nlReset(NL_TRUE) ;
633
634
635
	}
}

636
637
638
639
640
#ifndef DEBUG
Q_EXPORT_PLUGIN2(SurfaceDeformationPlugin, SurfaceDeformationPlugin)
#else
Q_EXPORT_PLUGIN2(SurfaceDeformationPluginD, SurfaceDeformationPlugin)
#endif
641
642
643
644

} // namespace SCHNApps

} // namespace CGoGN