map3MR_PrimalAdapt.hpp 31.1 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24
/*******************************************************************************
* CGoGN: Combinatorial and Geometric modeling with Generic N-dimensional Maps  *
* version 0.1                                                                  *
* Copyright (C) 2009-2012, IGG Team, LSIIT, University of Strasbourg           *
*                                                                              *
* This library is free software; you can redistribute it and/or modify it      *
* under the terms of the GNU Lesser General Public License as published by the *
* Free Software Foundation; either version 2.1 of the License, or (at your     *
* option) any later version.                                                   *
*                                                                              *
* This library is distributed in the hope that it will be useful, but WITHOUT  *
* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or        *
* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License  *
* for more details.                                                            *
*                                                                              *
* You should have received a copy of the GNU Lesser General Public License     *
* along with this library; if not, write to the Free Software Foundation,      *
* Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301 USA.           *
*                                                                              *
* Web site: http://cgogn.unistra.fr/                                           *
* Contact information: cgogn@unistra.fr                                        *
*                                                                              *
*******************************************************************************/

25
#include "Algo/Multiresolution/Map3MR/map3MR_PrimalAdapt.h"
26 27 28 29

namespace CGoGN
{

30
namespace Algo
31
{
32

33 34 35
namespace Volume
{

36 37 38 39 40 41 42 43 44 45 46 47
namespace MR
{

namespace Primal
{

namespace Adaptive
{

template <typename PFP>
Map3MR<PFP>::Map3MR(typename PFP::MAP& map) :
	m_map(map),
48
	shareVertexEmbeddings(false),
49 50 51 52 53 54
	vertexVertexFunctor(NULL),
	edgeVertexFunctor(NULL),
	faceVertexFunctor(NULL),
	volumeVertexFunctor(NULL)
{

55 56 57 58
}

/*! @name Cells informations
 *************************************************************************/
59 60
template <typename PFP>
unsigned int Map3MR<PFP>::edgeLevel(Dart d)
61
{
62
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"edgeLevel : called with a dart inserted after current level") ;
63

64 65
	// the level of an edge is the maximum of the
	// insertion levels of its darts
untereiner's avatar
untereiner committed
66 67 68 69 70 71 72
//	unsigned int r = 0;
//	Dart dit = d;
//	do
//	{
		unsigned int ld = m_map.getDartLevel(d) ;
		unsigned int ldd = m_map.getDartLevel(m_map.phi2(d)) ;
		unsigned int max = ld > ldd ? ld : ldd;
73

untereiner's avatar
untereiner committed
74
//		r =  r < max ? max : r ;
75

untereiner's avatar
untereiner committed
76 77
//		dit = m_map.alpha2(dit);
//	} while(dit != d);
78

untereiner's avatar
untereiner committed
79 80
//	return r;
	return max;
81 82
}

83 84
template <typename PFP>
unsigned int Map3MR<PFP>::faceLevel(Dart d)
85
{
86
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"faceLevel : called with a dart inserted after current level") ;
87

88
	if(m_map.getCurrentLevel() == 0)
89 90 91
		return 0 ;

	Dart it = d ;
92 93 94
	unsigned int min1 = m_map.getDartLevel(it) ;		// the level of a face is the second minimum of the
	it = m_map.phi1(it) ;
	unsigned int min2 = m_map.getDartLevel(it) ;		// insertion levels of its darts
95 96 97 98 99 100 101 102

	if(min2 < min1)
	{
		unsigned int tmp = min1 ;
		min1 = min2 ;
		min2 = tmp ;
	}

103
	it = m_map.phi1(it) ;
104 105
	while(it != d)
	{
106
		unsigned int dl = m_map.getDartLevel(it) ;
107 108 109 110 111 112 113 114 115 116
		if(dl < min2)
		{
			if(dl < min1)
			{
				min2 = min1 ;
				min1 = dl ;
			}
			else
				min2 = dl ;
		}
117
		it = m_map.phi1(it) ;
118 119 120 121 122
	}

	return min2 ;
}

123 124
template <typename PFP>
unsigned int Map3MR<PFP>::volumeLevel(Dart d)
125
{
126
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"volumeLevel : called with a dart inserted after current level") ;
127

128
	if(m_map.getCurrentLevel() == 0)
129 130 131 132 133
		return 0 ;

	unsigned int vLevel = std::numeric_limits<unsigned int>::max(); //hook sioux

	//First : the level of a volume is the minimum of the levels of its faces
134
	Traversor3WF<typename PFP::MAP> travF(m_map, d);
135 136 137 138 139 140 141 142 143 144 145 146 147 148
	for (Dart dit = travF.begin(); dit != travF.end(); dit = travF.next())
	{
		// in a first time, the level of a face
		//the level of the volume is the minimum of the
		//levels of its faces
		unsigned int fLevel = faceLevel(dit);
		vLevel = fLevel < vLevel ? fLevel : vLevel ;
	}

	//Second : the case of all faces regularly subdivided but not the volume itself

	return vLevel;
}

149 150
template <typename PFP>
Dart Map3MR<PFP>::faceOldestDart(Dart d)
151
{
152
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"faceOldestDart : called with a dart inserted after current level") ;
153 154 155

	Dart it = d ;
	Dart oldest = it ;
156
	unsigned int l_old = m_map.getDartLevel(oldest) ;
157 158
	do
	{
159
		unsigned int l = m_map.getDartLevel(it) ;
untereiner's avatar
untereiner committed
160 161 162
		if(l == 0)
			return it ;
		if(l < l_old)
163 164 165 166
		{
			oldest = it ;
			l_old = l ;
		}
167
		it = m_map.phi1(it) ;
168 169 170 171
	} while(it != d) ;
	return oldest ;
}

172 173
template <typename PFP>
Dart Map3MR<PFP>::volumeOldestDart(Dart d)
174
{
175
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"volumeOldestDart : called with a dart inserted after current level") ;
176 177

	Dart oldest = d;
178
	Traversor3WF<typename PFP::MAP> travF(m_map, d);
179 180 181 182
	for (Dart dit = travF.begin(); dit != travF.end(); dit = travF.next())
	{
		//for every dart in this face
		Dart old = faceOldestDart(dit);
183
		if(m_map.getDartLevel(old) < m_map.getDartLevel(oldest))
184 185 186 187 188 189
			oldest = old;
	}

	return oldest;
}

190 191
template <typename PFP>
bool Map3MR<PFP>::edgeIsSubdivided(Dart d)
192
{
193
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"edgeIsSubdivided : called with a dart inserted after current level") ;
194

195
	if(m_map.getCurrentLevel() == m_map.getMaxLevel())
196 197
		return false ;

198 199 200 201
	Dart d2 = m_map.phi2(d) ;
	m_map.incCurrentLevel() ;
	Dart d2_l = m_map.phi2(d) ;
	m_map.decCurrentLevel() ; ;
202 203 204 205 206 207
	if(d2 != d2_l)
		return true ;
	else
		return false ;
}

208 209
template <typename PFP>
bool Map3MR<PFP>::faceIsSubdivided(Dart d)
210
{
211
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"faceIsSubdivided : called with a dart inserted after current level") ;
212

213
	if(m_map.getCurrentLevel() == m_map.getMaxLevel())
214 215 216
		return false ;

	// a face whose level in the current level map is lower than
untereiner's avatar
untereiner committed
217
	// the current level can be subdivided to higher levels
218
	unsigned int fLevel = faceLevel(d) ;
219
	if(fLevel < m_map.getCurrentLevel())
220 221 222
		return false ;

	bool subd = false ;
untereiner's avatar
untereiner committed
223 224 225 226 227 228 229

	//Une face dont toute les aretes sont subdivise mais pas la face elle meme renvoie false .. sinon c'est true
	Dart dit = d;
	bool edgesAreSubdivided = true;
	do
	{
		edgesAreSubdivided &= edgeIsSubdivided(dit);
230
		dit = m_map.phi1(dit);
untereiner's avatar
untereiner committed
231 232 233 234
	}while(dit != d);

	if(edgesAreSubdivided)
	{
235
		m_map.incCurrentLevel() ;
236
		if(m_map.getDartLevel(m_map.phi1(m_map.phi1(d))) == m_map.getCurrentLevel()) //TODO a vérifier le phi1(phi1())
untereiner's avatar
untereiner committed
237
			subd = true ;
238
		m_map.decCurrentLevel() ;
untereiner's avatar
untereiner committed
239 240 241
	}
	else
		return false;
242 243 244 245

	return subd ;
}

246 247
template <typename PFP>
bool Map3MR<PFP>::volumeIsSubdivided(Dart d)
248
{
249
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"volumeIsSubdivided : called with a dart inserted after current level") ;
250 251

	unsigned int vLevel = volumeLevel(d);
252
	if(vLevel <= m_map.getCurrentLevel())
253 254 255 256
		return false;

	//Test if all faces are subdivided
	bool facesAreSubdivided = faceIsSubdivided(d) ;
257
	Traversor3WF<typename PFP::MAP> travF(m_map, d);
258 259 260 261 262 263 264
	for (Dart dit = travF.begin(); dit != travF.end(); dit = travF.next())
	{
		facesAreSubdivided &= faceIsSubdivided(dit) ;
	}

	//But not the volume itself
	bool subd = false;
265 266
	m_map.incCurrentLevel();
	if(facesAreSubdivided && m_map.getDartLevel(m_map.phi2(m_map.phi1(m_map.phi1(d)))) == m_map.getCurrentLevel())
267
		subd = true;
268
	m_map.decCurrentLevel() ;
269 270 271
	return subd;
}

272 273 274 275 276 277 278 279 280 281
/*! @name Topological helping functions
 *************************************************************************/
template <typename PFP>
void Map3MR<PFP>::swapEdges(Dart d, Dart e)
{
	if(!m_map.PFP::MAP::ParentMap::isBoundaryEdge(d) && !m_map.PFP::MAP::ParentMap::isBoundaryEdge(e))
	{
		Dart d2 = m_map.phi2(d);
		Dart e2 = m_map.phi2(e);

282
		m_map.PFP::MAP::ParentMap::swapEdges(d,e);
283

284 285 286 287 288
//		m_map.PFP::MAP::ParentMap::unsewFaces(d);
//		m_map.PFP::MAP::ParentMap::unsewFaces(e);
//
//		m_map.PFP::MAP::ParentMap::sewFaces(d, e);
//		m_map.PFP::MAP::ParentMap::sewFaces(d2, e2);
289 290 291 292 293 294 295 296 297 298 299 300 301 302

		if(m_map.template isOrbitEmbedded<VERTEX>())
		{
			m_map.template copyDartEmbedding<VERTEX>(d, m_map.phi2(m_map.phi_1(d)));
			m_map.template copyDartEmbedding<VERTEX>(e, m_map.phi2(m_map.phi_1(e)));
			m_map.template copyDartEmbedding<VERTEX>(d2, m_map.phi2(m_map.phi_1(d2)));
			m_map.template copyDartEmbedding<VERTEX>(e2, m_map.phi2(m_map.phi_1(e2)));
		}

		if(m_map.template isOrbitEmbedded<EDGE>())
		{

		}

303 304
//		if(m_map.template isOrbitEmbedded<VOLUME>())
//			m_map.template setOrbitEmbeddingOnNewCell<VOLUME>(d);
305 306


307 308 309 310
		m_map.duplicateDart(d);
		m_map.duplicateDart(d2);
		m_map.duplicateDart(e);
		m_map.duplicateDart(e2);
311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359
	}
}


template <typename PFP>
Dart Map3MR<PFP>::cutEdge(Dart d)
{
	Dart dit = d;
	do
	{
		Dart dd = m_map.phi2(dit) ;
		Dart d1 = m_map.phi1(dit);
		Dart dd1 = m_map.phi1(dd);

		m_map.duplicateDart(dit);
		m_map.duplicateDart(dd);
		m_map.duplicateDart(d1);
		m_map.duplicateDart(dd1);

		dit = m_map.alpha2(dit);
	}while(dit != d);

	Dart nd = m_map.cutEdge(d);

	return nd;
}

template <typename PFP>
void Map3MR<PFP>::splitFace(Dart d, Dart e)
{
	Dart dprev = m_map.phi_1(d) ;
	Dart eprev = m_map.phi_1(e) ;

	m_map.duplicateDart(d);
	m_map.duplicateDart(e);
	m_map.duplicateDart(dprev);
	m_map.duplicateDart(eprev);

	m_map.duplicateDart(m_map.phi3(d));
	m_map.duplicateDart(m_map.phi3(e));
	m_map.duplicateDart(m_map.phi3(dprev));
	m_map.duplicateDart(m_map.phi3(eprev));

	m_map.splitFace(d,e);
}

template <typename PFP>
void Map3MR<PFP>::splitVolume(std::vector<Dart>& vd)
{
360
	m_map.splitVolume(vd);
361 362 363 364

	for(std::vector<Dart>::iterator it = vd.begin() ; it != vd.end() ; ++it)
	{
		Dart dit = *it;
365
		m_map.duplicateDart(dit);
366

367 368
		Dart dit2 = m_map.phi2(dit);
		m_map.duplicateDart(dit2);
369

370 371
		Dart dit23 = m_map.phi3(dit2);
		m_map.duplicateDart(dit23);
372

373 374
		Dart dit232 = m_map.phi2(dit23);
		m_map.duplicateDart(dit232);
375 376 377 378
	}
}


379 380
/*! @name Subdivision
 *************************************************************************/
381 382
template <typename PFP>
void Map3MR<PFP>::subdivideEdge(Dart d)
383
{
384
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"subdivideEdge : called with a dart inserted after current level") ;
untereiner's avatar
untereiner committed
385
	assert(!edgeIsSubdivided(d) || !"Trying to subdivide an already subdivided edge") ;
386

387 388
	assert(m_map.getCurrentLevel() == edgeLevel(d) || !"Trying to subdivide an edge on a bad current level") ;

389
	m_map.incCurrentLevel();
untereiner's avatar
untereiner committed
390

391
	Dart nd = cutEdge(d);
392

393
	(*edgeVertexFunctor)(nd) ;
untereiner's avatar
untereiner committed
394

395
	m_map.decCurrentLevel() ;
396 397
}

untereiner's avatar
untereiner committed
398 399 400 401 402 403 404 405
template <typename PFP>
void Map3MR<PFP>::coarsenEdge(Dart d)
{
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"subdivideEdge : called with a dart inserted after current level") ;
	//assert(edgeCanBeCoarsened(d) || !"Trying to subdivide an already subdivided edge") ;

}

406
template <typename PFP>
407
void Map3MR<PFP>::subdivideFace(Dart d, bool triQuad)
408
{
409
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"subdivideFace : called with a dart inserted after current level") ;
410 411 412 413 414
	assert(!faceIsSubdivided(d) || !"Trying to subdivide an already subdivided face") ;

	unsigned int fLevel = faceLevel(d) ;
	Dart old = faceOldestDart(d) ;

415 416
	m_map.pushLevel() ;
	m_map.setCurrentLevel(fLevel) ;		// go to the level of the face to subdivide its edges
417 418 419 420 421 422 423 424 425

	unsigned int degree = 0 ;
	Dart it = old ;
	do
	{
		++degree ;						// compute the degree of the face

		if(!edgeIsSubdivided(it))
			subdivideEdge(it) ;			// and cut the edges (if they are not already)
426
		it = m_map.phi1(it) ;
427 428
	} while(it != old) ;

429
	m_map.setCurrentLevel(fLevel + 1) ;	// go to the next level to perform face subdivision
430

431 432 433 434 435 436 437 438 439 440 441 442 443 444 445 446 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 469 470 471 472 473 474
    if(triQuad && degree == 3)	// if subdividing a triangle
    {
        Dart dd = m_map.phi1(old) ;
        Dart e = m_map.phi1(dd) ;
        (*vertexVertexFunctor)(e) ;
        e = m_map.phi1(e) ;
        splitFace(dd, e) ;

        dd = e ;
        e = m_map.phi1(dd) ;
        (*vertexVertexFunctor)(e) ;
        e = m_map.phi1(e) ;
        splitFace(dd, e) ;

        dd = e ;
        e = m_map.phi1(dd) ;
        (*vertexVertexFunctor)(e) ;
        e = m_map.phi1(e) ;
        splitFace(dd, e) ;
    }
    else							// if subdividing a polygonal face
    {
        Dart dd = m_map.phi1(old) ;
        Dart next = m_map.phi1(dd) ;
        (*vertexVertexFunctor)(next) ;
        next = m_map.phi1(next) ;
        splitFace(dd, next) ;			// insert a first edge
        Dart ne = m_map.phi2(m_map.phi_1(dd));

        cutEdge(ne) ;					// cut the new edge to insert the central vertex

        dd = m_map.phi1(next) ;
        (*vertexVertexFunctor)(dd) ;
        dd = m_map.phi1(dd) ;
        while(dd != ne)					// turn around the face and insert new edges
        {								// linked to the central vertex
            splitFace(m_map.phi1(ne), dd) ;
            dd = m_map.phi1(dd) ;
            (*vertexVertexFunctor)(dd) ;
            dd = m_map.phi1(dd) ;
        }

        (*faceVertexFunctor)(m_map.phi1(ne)) ;
    }
475

476
	m_map.popLevel() ;
477 478
}

479
template <typename PFP>
480
unsigned int Map3MR<PFP>::subdivideVolume(Dart d, bool triQuad, bool OneLevelDifference)
481
{
482
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"subdivideVolume : called with a dart inserted after current level") ;
483 484 485 486 487
	assert(!volumeIsSubdivided(d) || !"Trying to subdivide an already subdivided face") ;

	unsigned int vLevel = volumeLevel(d);
	Dart old = volumeOldestDart(d);

488 489
	m_map.pushLevel() ;
	m_map.setCurrentLevel(vLevel) ;		// go to the level of the volume to subdivide its faces
490

491 492
	if(m_map.getCurrentLevel() == m_map.getMaxLevel())
		m_map.addLevelBack() ;
493 494

	//
untereiner's avatar
untereiner committed
495
	// Subdivide Faces and Edges
496
	//
497
	Traversor3WF<typename PFP::MAP> traF(m_map, old);
498 499 500 501
	for(Dart dit = traF.begin(); dit != traF.end(); dit = traF.next())
	{
		//if needed subdivide face
		if(!faceIsSubdivided(dit))
502
			subdivideFace(dit,triQuad);
503 504 505 506 507
	}

	//
	// Create inside volumes
	//
untereiner's avatar
untereiner committed
508 509
	std::vector<std::pair<Dart, Dart> > subdividedFaces;
	subdividedFaces.reserve(128);
510

511 512 513 514 515 516
	bool isocta = false;
	bool ishex = false;
	bool isprism = false;
	bool ispyra = false;

	Dart centralDart = NIL;
517
	Traversor3WV<typename PFP::MAP> traWV(m_map, d);
untereiner's avatar
untereiner committed
518
	for(Dart ditWV = traWV.begin(); ditWV != traWV.end(); ditWV = traWV.next())
519
	{
520
		m_map.incCurrentLevel() ;
521

untereiner's avatar
untereiner committed
522 523 524 525
		Dart e = ditWV;
		std::vector<Dart> v ;

		do
526
		{
527
			v.push_back(m_map.phi1(e));
528 529 530

			if(m_map.phi1(m_map.phi1(m_map.phi1(e))) != e)
				v.push_back(m_map.phi1(m_map.phi1(e)));
untereiner's avatar
untereiner committed
531

532
			if(!m_map.PFP::MAP::ParentMap::isBoundaryEdge(m_map.phi1(e)))
533
				subdividedFaces.push_back(std::pair<Dart,Dart>(m_map.phi1(e),m_map.phi2(m_map.phi1(e))));
untereiner's avatar
untereiner committed
534

535 536 537
			if(m_map.phi1(m_map.phi1(m_map.phi1(e))) != e)
				if(!m_map.PFP::MAP::ParentMap::isBoundaryEdge(m_map.phi1(m_map.phi1(e))))
					subdividedFaces.push_back(std::pair<Dart,Dart>(m_map.phi1(m_map.phi1(e)),m_map.phi2(m_map.phi1(m_map.phi1(e)))));
untereiner's avatar
untereiner committed
538

539
			e = m_map.phi2(m_map.phi_1(e));
540
		}
untereiner's avatar
untereiner committed
541
		while(e != ditWV);
542

543
		m_map.splitVolume(v);
544

545
		unsigned int fdeg = m_map.faceDegree(m_map.phi2(m_map.phi1(ditWV)));
546

547 548 549 550 551 552 553 554 555 556 557 558 559 560 561 562 563 564 565 566 567
		if(fdeg == 4)
		{
			if(m_map.PFP::MAP::ParentMap::vertexDegree(ditWV) == 3)
			{
				isocta = false;
				ispyra = true;

				Dart it = ditWV;
				if((m_map.faceDegree(it) == 3) && (m_map.faceDegree(m_map.phi2(it))) == 3)
				{
					it = m_map.phi2(m_map.phi_1(it));
				}
				else if((m_map.faceDegree(it) == 3) && (m_map.faceDegree(m_map.phi2(it)) == 4))
				{
					it = m_map.phi1(m_map.phi2(it));
				}

				Dart old = m_map.phi2(m_map.phi1(it));
				Dart dd = m_map.phi1(m_map.phi1(old));

				m_map.splitFace(old,dd) ;
568
				centralDart = old;
569 570 571
			}
			else
			{
572 573 574 575
				if(ispyra)
					isocta = false;
				else
					isocta = true;
576 577 578 579 580 581 582 583 584 585 586 587 588 589 590

				Dart old = m_map.phi2(m_map.phi1(ditWV));
				Dart dd = m_map.phi1(old) ;
				m_map.splitFace(old,dd) ;

				Dart ne = m_map.phi1(old);

				m_map.cutEdge(ne);
				centralDart = m_map.phi1(ne);

				Dart stop = m_map.phi2(m_map.phi1(ne));
				ne = m_map.phi2(ne);
				do
				{
					dd = m_map.phi1(m_map.phi1(ne));
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
					m_map.splitFace(ne, dd) ;

					ne = m_map.phi2(m_map.phi_1(ne));
					dd = m_map.phi1(dd);
				}
				while(dd != stop);
			}
		}
		else if(fdeg == 5)
		{
			isprism = true;

			Dart it = ditWV;
			if(m_map.faceDegree(it) == 3)
			{
				it = m_map.phi2(m_map.phi_1(it));
			}
			else if(m_map.faceDegree(m_map.phi2(m_map.phi_1(ditWV))) == 3)
			{
				it = m_map.phi2(m_map.phi_1(m_map.phi2(m_map.phi_1(it))));
			}

			Dart old = m_map.phi2(m_map.phi1(it));
			Dart dd = m_map.phi_1(m_map.phi_1(old));

			m_map.splitFace(old,dd) ;
		}
		if(fdeg == 6)
620
		{
621 622 623 624 625 626 627 628 629 630 631 632 633 634 635 636 637
			ishex = true;

			Dart dd = m_map.phi2(m_map.phi1(ditWV));;
			Dart next = m_map.phi1(m_map.phi1(dd)) ;
			m_map.splitFace(dd, next) ;		// insert a first edge

			Dart ne = m_map.phi2(m_map.phi_1(dd)) ;
			m_map.cutEdge(ne) ;				// cut the new edge to insert the central vertex
			centralDart = m_map.phi1(ne);

			dd = m_map.phi1(m_map.phi1(next)) ;
			while(dd != ne)				// turn around the face and insert new edges
			{							// linked to the central vertex
				Dart tmp = m_map.phi1(ne) ;
				m_map.splitFace(tmp, dd) ;
				dd = m_map.phi1(m_map.phi1(dd)) ;
			}
638
		}
untereiner's avatar
untereiner committed
639

640
		m_map.decCurrentLevel() ;
641 642
	}

643
	if(ishex)
644
	{
645 646 647 648 649 650 651 652 653 654 655 656 657 658 659 660
		m_map.incCurrentLevel();

		m_map.deleteVolume(m_map.phi3(m_map.phi2(m_map.phi1(d))));

		for (std::vector<std::pair<Dart,Dart> >::iterator it = subdividedFaces.begin(); it != subdividedFaces.end(); ++it)
		{
			Dart f1 = m_map.phi2((*it).first);
			Dart f2 = m_map.phi2((*it).second);

			if(m_map.isBoundaryFace(f1) && m_map.isBoundaryFace(f2))
			{
					m_map.sewVolumes(f1, f2);//, false);
			}
		}

		//replonger l'orbit de ditV.
661
		Algo::Topo::setOrbitEmbedding<VERTEX>(m_map,centralDart, m_map.template getEmbedding<VERTEX>(centralDart));
662 663 664
		(*volumeVertexFunctor)(centralDart) ;

		m_map.decCurrentLevel() ;
665
	}
666

667 668 669 670 671 672 673 674 675 676 677 678 679 680 681
	if(ispyra)
	{
		isocta = false;

		Dart ditV = d;

		Traversor3WV<typename PFP::MAP> traWV(m_map, d);
		for(Dart ditWV = traWV.begin(); ditWV != traWV.end(); ditWV = traWV.next())
		{
			if(m_map.PFP::MAP::ParentMap::vertexDegree(ditWV) == 4)
				ditV = ditWV;
		}

		m_map.incCurrentLevel();
		Dart x = m_map.phi_1(m_map.phi2(m_map.phi1(ditV)));
682
		std::vector<Dart> embVol;
683 684 685 686 687 688 689

		Dart f = x;
		do
		{
			Dart f3 = m_map.phi3(f);
			Dart tmp =  m_map.phi_1(m_map.phi2(m_map.phi_1(m_map.phi2(m_map.phi_1(f3))))); //future voisin par phi2
			swapEdges(f3, tmp);
690
			embVol.push_back(tmp);
691 692 693

			f = m_map.phi2(m_map.phi_1(f));
		}while(f != x);
untereiner's avatar
untereiner committed
694

695 696 697 698 699

		for(std::vector<Dart>::iterator it = embVol.begin() ; it != embVol.end() ; ++it)
		{
			Dart dit = *it;

700
			Algo::Topo::setOrbitEmbeddingOnNewCell<VOLUME>(m_map,dit);
701 702 703 704 705 706 707
			m_map.template copyCell<VOLUME>(dit, ditV);
		}

		//embed the volumes around swapEdges



708
		//replonger l'orbit de ditV.
709 710
		Algo::Topo::setOrbitEmbedding<VERTEX>(m_map, m_map.phi2(m_map.phi3(x)), m_map.template getEmbedding<VERTEX>(m_map.phi2(m_map.phi3(x))));

711 712
		//m_map.template setOrbitEmbedding<VERTEX>(centralDart, m_map.template getEmbedding<VERTEX>(centralDart));
		//(*volumeVertexFunctor)(x) ;
713 714 715 716 717

		m_map.decCurrentLevel() ;
	}

	if(isocta)
untereiner's avatar
untereiner committed
718
	{
719 720 721 722 723 724
		Traversor3WV<typename PFP::MAP> traWV(m_map, d);

		for(Dart ditWV = traWV.begin(); ditWV != traWV.end(); ditWV = traWV.next())
		{
			m_map.incCurrentLevel();
			Dart x = m_map.phi_1(m_map.phi2(m_map.phi1(ditWV)));
725
			std::vector<Dart> embVol;
726 727 728 729 730 731 732 733 734 735 736 737 738 739 740 741 742 743 744 745

			if(!Algo::Volume::Modelisation::Tetrahedralization::isTetrahedron<PFP>(m_map,x))
			{
				DartMarkerStore me(m_map);

				Dart f = x;

				do
				{
					Dart f3 = m_map.phi3(f);

					if(!me.isMarked(f3))
					{
						Dart tmp =  m_map.phi_1(m_map.phi2(m_map.phi_1(m_map.phi2(m_map.phi_1(f3))))); //future voisin par phi2

						Dart f32 = m_map.phi2(f3);
						swapEdges(f3, tmp);

						me.markOrbit<EDGE>(f3);
						me.markOrbit<EDGE>(f32);
746 747

						embVol.push_back(tmp);
748 749 750 751 752 753 754
					}

					f = m_map.phi2(m_map.phi_1(f));
				}while(f != x);

			}

755 756 757
			for(std::vector<Dart>::iterator it = embVol.begin() ; it != embVol.end() ; ++it)
			{
				Dart dit = *it;
758
				Algo::Topo::setOrbitEmbeddingOnNewCell<VOLUME>(m_map,dit);
759 760 761
				m_map.template copyCell<VOLUME>(dit, d);
			}

762
			Algo::Topo::setOrbitEmbeddingOnNewCell<VOLUME>(m_map, x, m_map.template getEmbedding<VERTEX>(x));
763
			(*volumeVertexFunctor)(x) ;
764 765 766

			m_map.decCurrentLevel() ;
		}
untereiner's avatar
untereiner committed
767 768
	}

769
	if(isprism)
untereiner's avatar
untereiner committed
770
	{
771 772 773 774 775 776 777 778 779 780 781 782 783 784 785 786 787 788 789 790 791 792 793 794 795 796 797 798 799 800 801 802 803 804 805 806 807 808 809
		m_map.incCurrentLevel();

		Dart ditWV = d;
		if(m_map.faceDegree(d) == 3)
		{
			ditWV = m_map.phi2(m_map.phi_1(d));
		}
		else if(m_map.faceDegree(m_map.phi2(m_map.phi_1(d))) == 3)
		{
			ditWV = m_map.phi1(m_map.phi2(d));
		}

		ditWV = m_map.phi3(m_map.phi_1(m_map.phi2(m_map.phi1(ditWV))));


		std::vector<Dart> path;

		Dart dtemp = ditWV;
		do
		{
			//future voisin par phi2
			Dart sF1 = m_map.phi1(m_map.phi2(m_map.phi3(dtemp)));
			Dart wrongVolume = m_map.phi3(sF1);
			Dart sF2 = m_map.phi3(m_map.phi2(wrongVolume));
			Dart tmp =  m_map.phi3(m_map.phi2(m_map.phi1(sF2)));
			swapEdges(dtemp, tmp);

			m_map.deleteVolume(wrongVolume);
			m_map.sewVolumes(sF1,sF2);

			path.push_back(dtemp);
			dtemp = m_map.phi_1(m_map.phi2(m_map.phi_1(dtemp)));


		}while(dtemp != ditWV);

		m_map.splitVolume(path);

		m_map.decCurrentLevel() ;
untereiner's avatar
untereiner committed
810 811
	}

812 813
	m_map.incCurrentLevel();
	m_map.popLevel() ;
814 815

	return vLevel;
816 817
}

818 819 820 821 822
template <typename PFP>
unsigned int Map3MR<PFP>::subdivideHexa(Dart d, bool OneLevelDifference)
{
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"subdivideVolume : called with a dart inserted after current level") ;
	assert(!volumeIsSubdivided(d) || !"Trying to subdivide an already subdivided face") ;
823

824 825
	unsigned int vLevel = volumeLevel(d);
	Dart old = volumeOldestDart(d);
826

827 828 829 830 831 832
	m_map.pushLevel() ;
	m_map.setCurrentLevel(vLevel) ;		// go to the level of the volume to subdivide its faces

	if(m_map.getCurrentLevel() == m_map.getMaxLevel())
		m_map.addLevelBack() ;

untereiner's avatar
untereiner committed
833
	if(OneLevelDifference)
834
	{
untereiner's avatar
untereiner committed
835 836
		Traversor3WF<typename PFP::MAP> traF(m_map, old);
		for(Dart dit = traF.begin(); dit != traF.end(); dit = traF.next())
837
		{
untereiner's avatar
untereiner committed
838 839 840 841
			Dart nv = m_map.phi3(dit);
			if(!m_map.isBoundaryMarked3(nv))
				if(volumeLevel(nv) == vLevel - 1)
					subdivideHexa(nv,false);
842 843
		}
	}
untereiner's avatar
untereiner committed
844 845
	std::vector<std::pair<Dart, Dart> > subdividedFaces;
	subdividedFaces.reserve(128);
846

847 848 849 850 851 852
	Traversor3WF<typename PFP::MAP> traF(m_map, old);
	for(Dart dit = traF.begin(); dit != traF.end(); dit = traF.next())
	{
		//if needed subdivide face
		if(!faceIsSubdivided(dit))
			subdivideFace(dit,false);
untereiner's avatar
untereiner committed
853 854 855 856 857 858 859 860 861 862 863 864 865 866 867 868 869

		//save a dart from the subdivided face
		m_map.pushLevel();

		unsigned int fLevel = faceLevel(dit); //puisque dans tous les cas, la face est subdivisee
		m_map.setCurrentLevel(fLevel + 1) ;

		//le brin est forcement du niveau cur
		Dart cf = m_map.phi1(dit);
		Dart e = cf;
		do
		{
			subdividedFaces.push_back(std::pair<Dart,Dart>(e,m_map.phi2(e)));
			e = m_map.phi2(m_map.phi1(e));
		}while (e != cf);

		m_map.popLevel();
870 871
	}

872

untereiner's avatar
untereiner committed
873
	Dart centralDart = NIL;
874

untereiner's avatar
untereiner committed
875
	Traversor3WV<typename PFP::MAP> traWV(m_map, old);
876 877
	for(Dart ditWV = traWV.begin(); ditWV != traWV.end(); ditWV = traWV.next())
	{
untereiner's avatar
untereiner committed
878 879 880
		m_map.pushLevel();

		m_map.setCurrentLevel(vLevel + 1) ;
881

882 883 884 885 886 887 888 889 890 891
		(*vertexVertexFunctor)(ditWV) ;

		Dart e = ditWV;
		std::vector<Dart> v ;

		do
		{
			v.push_back(m_map.phi1(e));
			v.push_back(m_map.phi1(m_map.phi1(e)));

untereiner's avatar
untereiner committed
892 893
//			if(!m_map.PFP::MAP::ParentMap::isBoundaryEdge(m_map.phi1(e)))
//				subdividedFaces.push_back(std::pair<Dart,Dart>(m_map.phi1(e),m_map.phi2(m_map.phi1(e))));
894

untereiner's avatar
untereiner committed
895 896 897
//			//if(m_map.phi1(m_map.phi1(m_map.phi1(e))) != e)
//			if(!m_map.PFP::MAP::ParentMap::isBoundaryEdge(m_map.phi1(m_map.phi1(e))))
//				subdividedFaces.push_back(std::pair<Dart,Dart>(m_map.phi1(m_map.phi1(e)),m_map.phi2(m_map.phi1(m_map.phi1(e)))));
898 899 900 901 902 903 904 905 906 907 908 909 910 911 912 913 914 915

			e = m_map.phi2(m_map.phi_1(e));
		}
		while(e != ditWV);

		m_map.splitVolume(v);

		Dart dd = m_map.phi2(m_map.phi1(ditWV));;
		Dart next = m_map.phi1(m_map.phi1(dd)) ;
		m_map.splitFace(dd, next) ;		// insert a first edge

		Dart ne = m_map.phi2(m_map.phi_1(dd)) ;
		m_map.cutEdge(ne) ;				// cut the new edge to insert the central vertex
		centralDart = m_map.phi1(ne);

		dd = m_map.phi1(m_map.phi1(next)) ;
		while(dd != ne)				// turn around the face and insert new edges
		{							// linked to the central vertex
untereiner's avatar
untereiner committed
916
			m_map.splitFace(m_map.phi1(ne), dd) ;
917 918 919
			dd = m_map.phi1(m_map.phi1(dd)) ;
		}

untereiner's avatar
untereiner committed
920
		m_map.popLevel() ;
921 922
	}

untereiner's avatar
untereiner committed
923
	m_map.setCurrentLevel(vLevel + 1) ;
924

untereiner's avatar
untereiner committed
925
	m_map.deleteVolume(m_map.phi3(m_map.phi2(m_map.phi1(old))));
926 927 928 929 930 931 932 933 934 935 936 937 938 939

	for (std::vector<std::pair<Dart,Dart> >::iterator it = subdividedFaces.begin(); it != subdividedFaces.end(); ++it)
	{
		Dart f1 = m_map.phi2((*it).first);
		Dart f2 = m_map.phi2((*it).second);

		if(m_map.isBoundaryFace(f1) && m_map.isBoundaryFace(f2))
		{
			m_map.sewVolumes(f1, f2);//, false);
		}
	}

	(*volumeVertexFunctor)(centralDart) ;

940 941 942 943
	m_map.popLevel() ;

	return vLevel;
}
944

945 946 947
template <typename PFP>
void Map3MR<PFP>::subdivideFace2(Dart d, bool triQuad)
{
untereiner's avatar
untereiner committed
948 949 950 951 952
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"subdivideFace : called with a dart inserted after current level") ;
	assert(!faceIsSubdivided(d) || !"Trying to subdivide an already subdivided face") ;

	unsigned int fLevel = faceLevel(d) ;
	Dart old = faceOldestDart(d) ;
953

untereiner's avatar
untereiner committed
954 955
	m_map.pushLevel() ;
	m_map.setCurrentLevel(fLevel) ;		// go to the level of the face to subdivide its edges
956

untereiner's avatar
untereiner committed
957
	std::cout << "face Level = " << fLevel << std::endl;
958

untereiner's avatar
untereiner committed
959 960 961 962 963 964 965 966 967 968 969 970
	unsigned int degree = 0 ;
	Dart it = old ;
	do
	{
		++degree ;						// compute the degree of the face

		std::cout << "edge Level = " << edgeLevel(it) << std::endl;
		if(!edgeIsSubdivided(it))
		{
			std::cout << "cutEdge" << std::endl;
			subdivideEdge(it) ;			// and cut the edges (if they are not already)
		}
971

untereiner's avatar
untereiner committed
972 973
		it = m_map.phi1(it) ;
	} while(it != old) ;
974

untereiner's avatar
untereiner committed
975
//	m_map.setCurrentLevel(fLevel + 1) ;	// go to the next level to perform face subdivision
976 977 978 979 980 981 982 983 984 985 986 987 988 989 990 991 992 993 994 995 996 997 998 999 1000 1001 1002 1003 1004 1005 1006 1007 1008 1009 1010 1011 1012 1013 1014 1015 1016 1017 1018 1019 1020 1021

//    if(triQuad && degree == 3)	// if subdividing a triangle
//    {
//        Dart dd = m_map.phi1(old) ;
//        Dart e = m_map.phi1(dd) ;
//        (*vertexVertexFunctor)(e) ;
//        e = m_map.phi1(e) ;
//        splitFace(dd, e) ;

//        dd = e ;
//        e = m_map.phi1(dd) ;
//        (*vertexVertexFunctor)(e) ;
//        e = m_map.phi1(e) ;
//        splitFace(dd, e) ;

//        dd = e ;
//        e = m_map.phi1(dd) ;
//        (*vertexVertexFunctor)(e) ;
//        e = m_map.phi1(e) ;
//        splitFace(dd, e) ;
//    }
//    else							// if subdividing a polygonal face
//    {
//        Dart dd = m_map.phi1(old) ;
//        Dart next = m_map.phi1(dd) ;
//        (*vertexVertexFunctor)(next) ;
//        next = m_map.phi1(next) ;
//        splitFace(dd, next) ;			// insert a first edge
//        Dart ne = m_map.phi2(m_map.phi_1(dd));

//        cutEdge(ne) ;					// cut the new edge to insert the central vertex

//        dd = m_map.phi1(next) ;
//        (*vertexVertexFunctor)(dd) ;
//        dd = m_map.phi1(dd) ;
//        while(dd != ne)					// turn around the face and insert new edges
//        {								// linked to the central vertex
//            splitFace(m_map.phi1(ne), dd) ;
//            dd = m_map.phi1(dd) ;
//            (*vertexVertexFunctor)(dd) ;
//            dd = m_map.phi1(dd) ;
//        }

//        (*faceVertexFunctor)(m_map.phi1(ne)) ;
//    }

untereiner's avatar
untereiner committed
1022
	m_map.popLevel() ;
1023
}
1024 1025 1026 1027

template <typename PFP>
unsigned int Map3MR<PFP>::subdivideHexa2(Dart d, bool OneLevelDifference)
{
untereiner's avatar
untereiner committed
1028 1029
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"subdivideVolume : called with a dart inserted after current level") ;
	assert(!volumeIsSubdivided(d) || !"Trying to subdivide an already subdivided face") ;
1030

untereiner's avatar
untereiner committed
1031 1032
	unsigned int vLevel = volumeLevel(d);
	Dart old = volumeOldestDart(d);
1033

untereiner's avatar
untereiner committed
1034 1035
	m_map.pushLevel() ;
	m_map.setCurrentLevel(vLevel) ;		// go to the level of the volume to subdivide its faces
1036

untereiner's avatar
untereiner committed
1037
	std::cout << "volume Level = " << vLevel << std::endl;
1038

untereiner's avatar
untereiner committed
1039 1040 1041 1042 1043
	if(m_map.getCurrentLevel() == m_map.getMaxLevel())
	{
		std::cout << "addLevel ?" << std::endl;
		m_map.addLevelBack() ;
	}
1044

untereiner's avatar
untereiner committed
1045 1046 1047 1048
	Traversor3WF<typename PFP::MAP> traF(m_map, old);
	for(Dart dit = traF.begin(); dit != traF.end(); dit = traF.next())
	{
		std::cout << "face" << std::endl;
1049

untereiner's avatar
untereiner committed
1050 1051 1052 1053 1054 1055 1056
		//if needed subdivide face
		if(!faceIsSubdivided(dit))
		{
			std::cout << "subdivide Face"<< std::endl;
			subdivideFace2(dit,false);
		}
	}
1057

untereiner's avatar
untereiner committed
1058
	m_map.popLevel() ;
1059

untereiner's avatar
untereiner committed
1060
	std::cout << "plop" << std::endl;
1061

untereiner's avatar
untereiner committed
1062
	return vLevel;
1063 1064 1065
}


1066 1067
template <typename PFP>
void Map3MR<PFP>::subdivideVolumeTetOcta(Dart d)
1068
{
1069
	assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"subdivideVolumeTetOcta : called with a dart inserted after current level") ;
1070 1071 1072 1073 1074
	assert(!volumeIsSubdivided(d) || !"Trying to subdivide an already subdivided face") ;

	unsigned int vLevel = volumeLevel(d);
	Dart old = volumeOldestDart(d);

1075 1076
	m_map.pushLevel() ;
	m_map.setCurrentLevel(vLevel) ;		// go to the level of the face to subdivide its edges
1077

1078 1079
	if(m_map.getCurrentLevel() == m_map.getMaxLevel())
		m_map.addLevelBack() ;
1080

untereiner's avatar
untereiner committed
1081 1082
	unsigned int j = 0;

1083
	//
untereiner's avatar
untereiner committed
1084
	// Subdivide Faces and Edges
1085
	//
1086
	Traversor3WF<typename PFP::MAP> traF(m_map, old);
1087 1088
	for(Dart dit = traF.begin(); dit != traF.end(); dit = traF.next())
	{
1089
		std::cout << "CurrentLevel = " << m_map.getCurrentLevel() << std::endl;
untereiner's avatar
untereiner committed
1090 1091 1092
		std::cout << "face level = " << faceLevel(dit) << std::endl;
		//std::cout << "d = " << dit << " is Subdivided ? " << faceIsSubdivided(dit) << std::endl;

1093 1094
		//if needed subdivide face
		if(!faceIsSubdivided(dit))
untereiner's avatar
untereiner committed
1095 1096
		{
			std::cout << "subdivide face = " << dit << std::endl;
1097
			subdivideFace(dit, true);
untereiner's avatar
untereiner committed
1098 1099
			++j;
		}
1100 1101 1102 1103 1104 1105 1106
	}

	//
	// Create inside volumes
	//
	Dart centralDart = NIL;
	bool isNotTet = false;
1107 1108
	Traversor3WV<typename PFP::MAP> traV(m_map, old);
	m_map.setCurrentLevel(vLevel + 1) ;
1109 1110
	for(Dart dit = traV.begin(); dit != traV.end(); dit = traV.next())
	{
1111
		//m_map.template setOrbitEmbedding<VERTEX>(dit, EMBNULL);
1112 1113
		(*vertexVertexFunctor)(dit) ;