surfaceDeformation.cpp 19.4 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
	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);
}

444
void SurfaceDeformationPlugin::changeSelectedMap(View* view, MapHandlerGen* map, bool fromUI)
445
446
447
448
449
450
451
452
453
{
	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
		if(!fromUI)
			m_dockTab->refreshUI(params);
463
464
465
	}
}

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

	if(view->isCurrentView())
474
475
476
477
	{
		if(!fromUI)
			m_dockTab->refreshUI(params);
	}
478
479
}

480
void SurfaceDeformationPlugin::changeVerticesSelectionMode(View* view, MapHandlerGen* map, SelectionMode m, bool fromUI)
481
482
483
484
485
{
	ParameterSet* params = h_viewParams[view];
	params->perMap[map->getName()]->verticesSelectionMode = m;

	if(view->isCurrentView())
486
487
488
489
	{
		if(!fromUI)
			m_dockTab->refreshUI(params);
	}
490
491
}

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

497
498
499
	nlMakeCurrent(perMap->nlContext);
	if(nlGetCurrentState() == NL_STATE_INITIAL)
		nlBegin(NL_SYSTEM) ;
500
501
	for(int coord = 0; coord < 3; ++coord)
	{
502
503
		LinearSolving::setupVariables<PFP2>(*map, perMap->vIndex, *perMap->lockingMarker, perMap->positionAttribute, coord) ;
		nlBegin(NL_MATRIX) ;
504
505
		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) ;
506
507
508
509
510
		nlEnd(NL_MATRIX) ;
		nlEnd(NL_SYSTEM) ;
		nlSolve() ;
		LinearSolving::getResult<PFP2>(*map, perMap->vIndex, perMap->positionAttribute, coord) ;
		nlReset(NL_TRUE) ;
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
	}
}

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] ;
530
//			PFP2::REAL area = perMap->vertexArea[d] ;
531
532
533
534
535
536
537
538
			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)
539
						cov(i,j) += v[i] * vv[j];// * perMap->edgeWeight[it] / area ;
540
541
542
543
544
545
546
				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)
547
							cov(i,j) += v[i] * vv[j];// * perMap->edgeWeight[dboundary] / area ;
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
619
620
621
622
623
624
625
				}
				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]) ;
		}
	}

626
627
628
	nlMakeCurrent(perMap->nlContext);
	if(nlGetCurrentState() == NL_STATE_INITIAL)
		nlBegin(NL_SYSTEM) ;
629
630
	for(int coord = 0; coord < 3; ++coord)
	{
631
632
		LinearSolving::setupVariables<PFP2>(*map, perMap->vIndex, *perMap->lockingMarker, perMap->positionAttribute, coord) ;
		nlBegin(NL_MATRIX) ;
633
634
//		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) ;
635
636
637
638
639
		nlEnd(NL_MATRIX) ;
		nlEnd(NL_SYSTEM) ;
		nlSolve() ;
		LinearSolving::getResult<PFP2>(*map, perMap->vIndex, perMap->positionAttribute, coord) ;
		nlReset(NL_TRUE) ;
640
641
642
	}
}

643
644
645
646
647
#ifndef DEBUG
Q_EXPORT_PLUGIN2(SurfaceDeformationPlugin, SurfaceDeformationPlugin)
#else
Q_EXPORT_PLUGIN2(SurfaceDeformationPluginD, SurfaceDeformationPlugin)
#endif
648
649
650
651

} // namespace SCHNApps

} // namespace CGoGN