-using namespace std;
+// MEFISTO2: a library to compute 2D triangulation from segmented boundaries
+//
+// Copyright (C) 2006 Laboratoire J.-L. Lions UPMC Paris
+//
+// 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.
+//
+// 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., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+// See http://www.ann.jussieu.fr/~perronne or email Perronnet@ann.jussieu.fr
+//
+// File : aptrte.cxx le C++ de l'appel du trianguleur plan
+// Module : SMESH
+// Author : Alain PERRONNET
+// Date : 16 mars 2006
+
#include "Rn.h"
#include "aptrte.h"
#include "utilities.h"
+using namespace std;
+
extern "C"
{
R aretemaxface_;
}
-void aptrte( Z nutysu, R aretmx,
- Z nblf, Z * nudslf, R2 * uvslf,
- Z nbpti, R2 *uvpti,
- Z & nbst, R2 * & uvst, Z & nbt, Z * & nust,
+void aptrte( Z nutysu, R aretmx,
+ Z nblf, Z * nudslf, R2 * uvslf,
+ Z nbpti, R2 * uvpti,
+ Z & nbst, R2 * & uvst,
+ Z & nbt, Z * & nust,
Z & ierr )
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// but : appel de la triangulation par un arbre-4 recouvrant
// ierr : 0 si pas d'erreur
// > 0 sinon
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-// auteur : Alain Perronnet Analyse Numerique Paris UPMC decembre 2001
+// auteur : Alain Perronnet Laboratoire J.-L. LIONS Paris UPMC mars 2006
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
{
+ Z nbsttria=4; //Attention: 4 sommets stockes par triangle
+ //no st1, st2, st3, 0 (non quadrangle)
+
R d, tcpu=0;
R3 direction=R3(0,0,0); //direction pour areteideale() inactive ici!
Z nbarfr=nudslf[nblf]; //nombre total d'aretes des lignes fermees
R3 *mnpxyd=NULL;
Z *mnsoar=NULL, mosoar=7, mxsoar, n1soar; //le hachage des aretes
Z *mnartr=NULL, moartr=3, mxartr, n1artr; //le no des 3 aretes des triangles
- Z *mntree=NULL, motree=9, mxtree; //L'arbre 4 de TE et nombre d'entiers par TE
+ Z *mntree=NULL, motree=9, mxtree; //L'arbre 4 de TE et nombre d'entiers par TE
Z *mnqueu=NULL, mxqueu;
Z *mn1arcf=NULL;
Z *mnarcf=NULL, mxarcf;
Z *mnarcf2=NULL;
Z *mnarcf3=NULL;
Z *mntrsu=NULL;
- Z *mndalf=NULL;
Z *mnslig=NULL;
Z *mnarst=NULL;
Z *mnlftr=NULL;
Z i, l, n, ns, ns0, ns1, ns2, nosotr[3], nt;
Z mxsomm, nbsomm, nbarpi, nbarli, ndtri0, mn;
Z moins1=-1;
+ R dist;
aretemaxface_ = aretmx;
// quelques reservations de tableaux pour faire les calculs
// ========================================================
- // le tableau pointeur sur la premiere arete de chaque ligne fermee
- if( mndalf!=NULL ) delete [] mndalf;
- mndalf = new Z[1+nblf];
- if( mndalf==NULL ) goto ERREUR;
- mndalf[0]=0;
-
// declaration du tableau des coordonnees des sommets de la frontiere
// puis des sommets internes ajoutes
// majoration empirique du nombre de sommets de la triangulation
- i = 4*nbarfr/10;
+ i = 4*nbarfr/10;
mxsomm = Max( 20000, 64*nbpti+i*i );
- MESSAGE( "APTRTE: Depart de la triangulation avec " );
+ MESSAGE( "APTRTE: Debut de la triangulation plane avec " );
MESSAGE( "nutysu=" << nutysu << " aretmx=" << aretmx << " mxsomm=" << mxsomm );
+ MESSAGE( nbarfr << " sommets sur la frontiere et " << nbpti << " points internes");
NEWDEPART:
//mnpxyd( 3, mxsomm ) les coordonnees UV des sommets et la taille d'arete aux sommets
mnpxyd[ns0].x = uvslf[ns0].x;
mnpxyd[ns0].y = uvslf[ns0].y;
mnpxyd[ns0].z = areteideale_( mnpxyd[ns0], direction );
-// cout << "Sommet " << ns0 << ": " << mnpxyd[ns0].x
-// << " " << mnpxyd[ns0].y << " longueur arete=" << mnpxyd[ns0].z << endl;
+// MESSAGE("Sommet " << ns0 << ": " << mnpxyd[ns0].x
+// << " " << mnpxyd[ns0].y << " longueur arete=" << mnpxyd[ns0].z);
//carre de la longueur de l'arete 1 de la ligne fermee n
d = pow( uvslf[ns0+1].x - uvslf[ns0].x, 2 )
//pas de test sur ierr car pas de saturation possible a ce niveau
//le pointeur dans le hachage sur la premiere arete de la ligne fermee n
- mndalf[n] = noar0;
+ //mndalf[n] = noar0;
//la nouvelle arete est la suivante de l'arete definie juste avant
if( noar > 0 )
mnpxyd[ns].x = uvslf[ns].x;
mnpxyd[ns].y = uvslf[ns].y;
mnpxyd[ns].z = areteideale_( mnpxyd[ns], direction );
-// cout << "Sommet " << ns << ": " << mnpxyd[ns].x
-// << " " << mnpxyd[ns].y << " longueur arete=" << mnpxyd[ns].z << endl;
+// MESSAGE("Sommet " << ns << ": " << mnpxyd[ns].x
+// << " " << mnpxyd[ns].y << " longueur arete=" << mnpxyd[ns].z);
//carre de la longueur de l'arete
d = pow( uvslf[ns2-1].x - uvslf[ns1-1].x, 2)
aremax = sqrt( aremax ); //longueur maximale d'une arete
aretmx = Min( aretmx, aremax ); //pour homogeneiser
- cout << "nutysu=" << nutysu << " aretmx=" << aretmx
- << " arete min=" << aremin << " arete max=" << aremax << endl;
+ MESSAGE("nutysu=" << nutysu << " aretmx=" << aretmx
+ << " arete min=" << aremin << " arete max=" << aremax);
//chainage des aretes frontalieres : la derniere arete frontaliere
mnsoar[ mosoar * noar - mosoar + 5 ] = 0;
deltacpu_( d );
tcpu += d;
- cout << "Temps de l'adaptation et l'homogeneisation de l'arbre-4 des TE="
- << d << " secondes" << endl;
+ MESSAGE("Temps de l'adaptation et l'homogeneisation de l'arbre-4 des TE="
+ << d << " secondes");
if( ierr != 0 )
{
//destruction du tableau auxiliaire et de l'arbre
MESSAGE( "Nombre d'echanges des diagonales de 2 triangles=" << n );
deltacpu_( d );
tcpu += d;
- cout << "Temps de la triangulation Delaunay par echange des diagonales="
- << d << " secondes" << endl;
+ MESSAGE("Temps de la triangulation Delaunay par echange des diagonales="
+ << d << " secondes");
//qualites de la triangulation actuelle
qualitetrte( mnpxyd, mosoar, mxsoar, mnsoar, moartr, mxartr, mnartr,
MESSAGE( "Restauration de " << n << " aretes perdues de la frontiere" );
deltacpu_( d );
tcpu += d;
- cout << "Temps de la recuperation des aretes perdues de la frontiere="
- << d << " secondes" << endl;
+ MESSAGE("Temps de la recuperation des aretes perdues de la frontiere="
+ << d << " secondes");
if( ierr != 0 ) goto ERREUR;
// mise en delaunay de la triangulation
// =====================================================
mnarcf3 = new Z[mxarcf];
- if( mnarcf3 == NULL ) goto ERREUR;
-
+ if( mnarcf3 == NULL )
+ {
+ cout << "aptrte: MC saturee mnarcf3=" << mnarcf3 << endl;
+ goto ERREUR;
+ }
teamqt_( nutysu,
mnarst, mosoar, mxsoar, n1soar, mnsoar,
moartr, mxartr, n1artr, mnartr,
mn1arcf, mnarcf, mnarcf1,
comxmi, nbarpi, nbsomm, mxsomm, mnpxyd, mnslig,
ierr );
+ if( mnarcf3 != NULL ) {delete [] mnarcf3; mnarcf3=NULL;}
if( mn1arcf != NULL ) {delete [] mn1arcf; mn1arcf=NULL;}
if( mnarcf != NULL ) {delete [] mnarcf; mnarcf =NULL;}
if( mnarcf1 != NULL ) {delete [] mnarcf1; mnarcf1=NULL;}
if( mnarcf2 != NULL ) {delete [] mnarcf2; mnarcf2=NULL;}
- if( mnarcf3 != NULL ) {delete [] mnarcf3; mnarcf3=NULL;}
deltacpu_( d );
tcpu += d;
// -----------------------------------------------------
// boucle sur les triangles occupes (internes et externes)
if( nust != NULL ) delete [] nust;
- nust = new Z[4*nbt];
+ nust = new Z[nbsttria*nbt];
if( nust == NULL ) goto ERREUR;
nbt = 0;
for (i=1; i<=mxartr; i++)
nust[nbt++] = 0;
}
}
- nbt /= 4; //le nombre final de triangles de la surface
- cout << "Nombre de sommets=" << nbst
- << " Nombre de triangles=" << nbt << endl;
-
+ nbt /= nbsttria; //le nombre final de triangles de la surface
+ MESSAGE( "APTRTE: Fin de la triangulation plane avec "<<nbst<<" sommets et "
+ << nbt << " triangles=" << nbt);
deltacpu_( d );
tcpu += d;
- MESSAGE( "Temps total de la triangulation=" << tcpu << " secondes" );
+ MESSAGE( "APTRTE: Temps total de la triangulation plane=" << tcpu << " secondes" );
// destruction des tableaux auxiliaires
// ------------------------------------
if( mnslig != NULL ) delete [] mnslig;
if( mnsoar != NULL ) delete [] mnsoar;
if( mnpxyd != NULL ) delete [] mnpxyd;
- if( mndalf != NULL ) delete [] mndalf;
if( mntree != NULL ) delete [] mntree;
if( mnqueu != NULL ) delete [] mnqueu;
if( mntrsu != NULL ) delete [] mntrsu;
}
else
{
- MESSAGE( "Triangulation non realisee " << ierr );
+ MESSAGE( "APTRTE: Triangulation NON REALISEE avec erreur=" << ierr );
if( ierr == 0 ) ierr=1;
goto NETTOYAGE;
}
{
//un triangle d'aire negative de plus
nbtrianeg++;
- cout << "ATTENTION: le triangle " << nt << " de sommets:"
+ MESSAGE("ATTENTION: le triangle " << nt << " de sommets:"
<< nosotr[0] << " " << nosotr[1] << " " << nosotr[2]
- << " a une aire " << d <<"<=0" << endl;
+ << " a une aire " << d <<"<=0");
}
//aire des triangles actuels
//les affichages
quamoy /= nbtria;
- cout << "Qualite moyenne=" << quamoy
+ MESSAGE("Qualite moyenne=" << quamoy
<< " Qualite minimale=" << quamin
- << " des " << nbtria << " triangles de surface totale="
- << aire << endl;
+ << " des " << nbtria << " triangles de surface plane totale="
+ << aire);
if( nbtrianeg>0 )
MESSAGE( "ATTENTION: nombre de triangles d'aire negative=" << nbtrianeg );