fct_generateur_frontiere.cpp

Aller à la documentation de ce fichier.

//####//------------------------------------------------------------
//####//------------------------------------------------------------
//####//  MAGiC
//####//  Jean Christophe Cuilliere et Vincent FRANCOIS
//####//  Departement de Genie Mecanique - UQTR
//####//------------------------------------------------------------
//####//  MAGIC est un projet de recherche de l equipe ERICCA
//####//  du departement de genie mecanique de l Universite du Quebec a Trois Rivieres
//####//  http://www.uqtr.ca/ericca
//####//  http://www.uqtr.ca/
//####//------------------------------------------------------------
//####//------------------------------------------------------------
//####//
//####//       fct_generateur_frontiere.cpp
//####//
//####//------------------------------------------------------------
//####//------------------------------------------------------------
//####//    COPYRIGHT 2000-2024
//####//  jeu 13 jun 2024 11:58:52 EDT
//####//------------------------------------------------------------
//####//------------------------------------------------------------
 
 
#include "gestionversion.h"
 
#pragma hdrstop
 
 
#include <math.h>
#include "fct_generateur_frontiere.h"
#include "mg_gestionnaire.h"
#include "ct_face.h"
#include "ct_point.h"
 
 
 
 
FCT_GENERATEUR_FRONTIERE::FCT_GENERATEUR_FRONTIERE(MG_GESTIONNAIRE& mggest,double dn,double eps,int ind1,double ind2,double *lst_pt_ctrl,int nbpt):FCT_GENERATEUR_ECHANTILLONAGE_ESTIMATION(dn,ind1,ind2),gest(mggest),epsilon(eps),point_controle(lst_pt_ctrl),nb_pt_ctrl(nbpt)
{
}
 
FCT_GENERATEUR_FRONTIERE::~FCT_GENERATEUR_FRONTIERE()
{
}
 
 
void FCT_GENERATEUR_FRONTIERE::echantillonner(MG_GROUPE_TOPOLOGIQUE* mggt)
{
    TPL_MAP_ENTITE<MG_ELEMENT_TOPOLOGIQUE*> lst;
    if (mggt!=NULL)
    {
      std::map<MG_ELEMENT_TOPOLOGIQUE*,MG_ELEMENT_TOPOLOGIQUE*>::iterator it;
      for(MG_ELEMENT_TOPOLOGIQUE* ele=mggt->get_premier(it);ele!=NULL;ele=mggt->get_suivant(it))
      {
        lst.ajouter(ele);
        ele->get_topologie_sousjacente(&lst);
      }
    }
    xmin=1e308,ymin=1e308,zmin=1e308;
    xmax=-1e308,ymax=-1e308,zmax=-1e308;
    int nb_geo=gest.get_nb_mg_geometrie();
    for (int i=0;i<nb_geo;i++)
    {
        MG_GEOMETRIE* mggeo=gest.get_mg_geometrie(i);
        int nb_face=mggeo->get_nb_mg_face();
        TPL_MAP_ENTITE<CT_FACE*> lstface;
        for (int j=0;j<nb_face;j++)
        {
            MG_FACE* face=mggeo->get_mg_face(j);
            CT_FACE* cface=new CT_FACE(face);
            cface->change_umax(-1e308);
            cface->change_umin(1e308);
            cface->change_vmax(-1e308);
            cface->change_vmin(1e308);
            lstface.ajouter(cface);
        }
        int nb_arete=mggeo->get_nb_mg_arete();
        for (int j=0;j<nb_arete;j++)
        {
            MG_ARETE* arete=mggeo->get_mg_arete(j);
            if (mggt!=NULL)
                if (lst.existe(arete)==0) continue;
            int nbfac=arete->get_nb_mg_coarete();
            double tmin=arete->get_tmin();
            double tmax=arete->get_tmax();
            double tdemi=0.5*(tmax+tmin);
            double xyz1[3],xyz2[3];
            arete->evaluer(tmin,xyz1);
            arete->evaluer(tdemi,xyz2);
            OT_VECTEUR_3D vec(xyz1,xyz2);
            double densdep=vec.get_longueur()*2.;
            MG_FACE** lstfac=new MG_FACE*[nbfac];
            for (int k=0;k<nbfac;k++)
                lstfac[k]=arete->get_mg_coarete(k)->get_boucle()->get_mg_face();
            if (pas_echantillon!=-1)
                for (int k=0;k<pas_echantillon+1;k++)
                {
                    double t=tmin+1.0*k/pas_echantillon*(tmax-tmin);
                    double xyz[3];
                    arete->evaluer(t,xyz);
                    if (xyz[0]<xmin) xmin=xyz[0];
                    if (xyz[0]>xmax) xmax=xyz[0];
                    if (xyz[1]<ymin) ymin=xyz[1];
                    if (xyz[1]>ymax) ymax=xyz[1];
                    if (xyz[2]<zmin) zmin=xyz[2];
                    if (xyz[2]>zmax) zmax=xyz[2];
                    double courbure=arete->get_M(t);
                    for (int l=0;l<nbfac;l++)
                    {
                        double uv[2];
                        lstfac[l]->inverser(uv,xyz);
                        CT_FACE* fac=lstface.getid(lstfac[l]->get_id());
                        if (uv[0]>fac->get_umax()) fac->change_umax(uv[0]);
                        if (uv[0]<fac->get_umin()) fac->change_umin(uv[0]);
                        if (uv[1]>fac->get_vmax()) fac->change_vmax(uv[1]);
                        if (uv[1]<fac->get_vmin()) fac->change_vmin(uv[1]);
                        double cmin,cmax;
                        lstfac[l]->get_courbure(uv,cmax,cmin);
                        cmax=fabs(cmax);
                        if (cmax>courbure) courbure=cmax;
                    }
                    double dn=densdep;
                    if (!(OPERATEUR::egal(courbure,0.0,1e-8))) dn=sqrt(8*epsilon/courbure);
                    else dn=dng;
                    if (dn>dng) dn=dng;
                    if (densdep<dn) dn=densdep;
                    double dist=zonezoom*(dng-dn);
                    if (!(OPERATEUR::egal(dist,0.0,1e-8)))
                    {
                        BOITE_3D boite(xyz[0]-dist,xyz[1]-dist,xyz[2]-dist,xyz[0]+dist,xyz[1]+dist,xyz[2]+dist);
                        DOUBLEN<1> ecart(&dn);
                        CT_POINT<1>* pt=new CT_POINT<1>(xyz[0],xyz[1],xyz[2],ecart,boite);
                        lstpoint.ajouter(pt);
                    }
                }
            else
            {
                double t=tmin;
                int ok=0;
                do
                {
                    double xyz[3];
                    arete->evaluer(t,xyz);
                    if (xyz[0]<xmin) xmin=xyz[0];
                    if (xyz[0]>xmax) xmax=xyz[0];
                    if (xyz[1]<ymin) ymin=xyz[1];
                    if (xyz[1]>ymax) ymax=xyz[1];
                    if (xyz[2]<zmin) zmin=xyz[2];
                    if (xyz[2]>zmax) zmax=xyz[2];
                    double courbure=arete->get_M(t);
                    for (int l=0;l<nbfac;l++)
                    {
                        double uv[2];
                        lstfac[l]->inverser(uv,xyz);
                        CT_FACE* fac=lstface.getid(lstfac[l]->get_id());
                        if (uv[0]>fac->get_umax()) fac->change_umax(uv[0]);
                        if (uv[0]<fac->get_umin()) fac->change_umin(uv[0]);
                        if (uv[1]>fac->get_vmax()) fac->change_vmax(uv[1]);
                        if (uv[1]<fac->get_vmin()) fac->change_vmin(uv[1]);
                        double cmin,cmax;
                        lstfac[l]->get_courbure(uv,cmax,cmin);
                        cmax=fabs(cmax);
                        if (cmax>courbure) courbure=cmax;
                    }
                    double dn=densdep;
                    double dntmp;
                    if (!(OPERATEUR::egal(courbure,0.0,1e-8))) dntmp=sqrt(8*epsilon/courbure);
                    else dntmp=1e308;
                    if (dntmp>dng) dn=dng;
                    else dn=dntmp;
                    if (densdep<dn) dn=densdep;
                    double dist=zonezoom*(dng-dn);
                    if (!(OPERATEUR::egal(dist,0.0,1e-8)))
                    {
                        BOITE_3D boite(xyz[0]-dist,xyz[1]-dist,xyz[2]-dist,xyz[0]+dist,xyz[1]+dist,xyz[2]+dist);
                        DOUBLEN<1> ecart(&dn);
                        CT_POINT<1>* pt=new CT_POINT<1>(xyz[0],xyz[1],xyz[2],ecart,boite);
                        lstpoint.ajouter(pt);
                    }
                    double dxyz[3];
                    arete->deriver(t,dxyz);
                    t=t+dntmp/sqrt(dxyz[0]*dxyz[0]+dxyz[1]*dxyz[1]+dxyz[2]*dxyz[2]);
                    if (t>tmax)
                    {
                        t=tmax;
                        ok++;
                    }
                }
                while (ok<2);
            }
            delete [] lstfac;
        }
 
        for (int j=0;j<nb_face;j++)
        {
            MG_FACE* face=mggeo->get_mg_face(j);
            if (mggt!=NULL)
                if (lst.existe(face)==0) continue;
            CT_FACE* cface=lstface.getid(face->get_id());
            double umin=cface->get_umin();
            double vmin=cface->get_vmin();
            double umax=cface->get_umax();
            double vmax=cface->get_vmax();
            if (pas_echantillon!=-1)
            {
                for (int k=0;k<pas_echantillon+1;k++)
                    for (int l=0;l<pas_echantillon+1;l++)
                    {
                        double uv[2];
                        uv[0]=umin+1.0*k/pas_echantillon*(umax-umin);
                        uv[1]=vmin+1.0*l/pas_echantillon*(vmax-vmin);
                        double xyz[3];
                        face->evaluer(uv,xyz);
                        if (xyz[0]<xmin) xmin=xyz[0];
                        if (xyz[1]<ymin) ymin=xyz[1];
                        if (xyz[2]<zmin) zmin=xyz[2];
                        if (xyz[0]>xmax) xmax=xyz[0];
                        if (xyz[1]>ymax) ymax=xyz[1];
                        if (xyz[2]>zmax) zmax=xyz[2];
                        double cmax,cmin;
                        face->get_courbure(uv,cmax,cmin);
                        cmax=fabs(cmax);
                        double dn;
                        if (!(OPERATEUR::egal(cmax,0.0,1e-8))) dn=sqrt(8*epsilon/cmax);
                        else dn=dng;
                        if (dn>dng) dn=dng;
                        double dist=zonezoom*(dng-dn);
                        if (!(OPERATEUR::egal(dist,0.0,1e-8)))
                        {
                            BOITE_3D boite(xyz[0]-dist,xyz[1]-dist,xyz[2]-dist,xyz[0]+dist,xyz[1]+dist,xyz[2]+dist);
                            DOUBLEN<1> ecart(&dn);
                            CT_POINT<1>* pt=new CT_POINT<1>(xyz[0],xyz[1],xyz[2],ecart,boite);
                            lstpoint.ajouter(pt);
                        }
                    }
            }
            else
            {
                double uv[2];
                double udepnou;
                double vdepnou;
                double udep=umin;
                double vdep=vmin;
                uv[0]=udep;
                uv[1]=vdep;
                double xyz[3];
                int ok=0;
                do
                {
                    face->evaluer(uv,xyz);
                    if (xyz[0]<xmin) xmin=xyz[0];
                    if (xyz[1]<ymin) ymin=xyz[1];
                    if (xyz[2]<zmin) zmin=xyz[2];
                    if (xyz[0]>xmax) xmax=xyz[0];
                    if (xyz[1]>ymax) ymax=xyz[1];
                    if (xyz[2]>zmax) zmax=xyz[2];
                    double cmax,cmin;
                    face->get_courbure(uv,cmax,cmin);
                    cmax=fabs(cmax);
                    double dn,dntmp;
                    if (!(OPERATEUR::egal(cmax,0.0,1e-8))) dntmp=sqrt(8*epsilon/cmax);
                    else dntmp=1e308;
                    if (dntmp>dng) dn=dng;
                    else dn=dntmp;
                    double dist=zonezoom*(dng-dn);
                    if (!(OPERATEUR::egal(dist,0.0,1e-8)))
                    {
                        BOITE_3D boite(xyz[0]-dist,xyz[1]-dist,xyz[2]-dist,xyz[0]+dist,xyz[1]+dist,xyz[2]+dist);
                        DOUBLEN<1> ecart(&dn);
                        CT_POINT<1>* pt=new CT_POINT<1>(xyz[0],xyz[1],xyz[2],ecart,boite);
                        lstpoint.ajouter(pt);
                    }
                    double dxyzdu[3];
                    double dxyzdv[3];
                    face->deriver(uv,dxyzdu,dxyzdv);
                    if (ok==-1) ok=6;
                    if (ok==0)
                    {
                        udepnou=uv[0]+dntmp/sqrt(dxyzdu[0]*dxyzdu[0]+dxyzdu[1]*dxyzdu[1]+dxyzdu[2]*dxyzdu[2]);
                        vdepnou=uv[1]+dntmp/sqrt(dxyzdv[0]*dxyzdv[0]+dxyzdv[1]*dxyzdv[1]+dxyzdv[2]*dxyzdv[2]);
                        if (udepnou>umax) udepnou=umax;
                        if (vdepnou>vmax) vdepnou=vmax;
                        ok=1;
                    }
                    if (ok==1)
                    {
                        uv[0]=uv[0]+dntmp/sqrt(dxyzdu[0]*dxyzdu[0]+dxyzdu[1]*dxyzdu[1]+dxyzdu[2]*dxyzdu[2]);
                        if (uv[0]>umax)
                        {
                            uv[0]=umax;
                            ok=2;
                        }
                    }
                    else if (ok==2)
                    {
                        uv[0]=udep;
                        uv[1]=vdepnou;
                        ok=3;
                        if (uv[1]+1e-10>vmax) ok=4;
                    }
                    else if (ok==3)
                    {
                        uv[1]=uv[1]+dntmp/sqrt(dxyzdv[0]*dxyzdv[0]+dxyzdv[1]*dxyzdv[1]+dxyzdv[2]*dxyzdv[2]);
                        if (uv[1]>vmax)
                        {
                            uv[1]=vmax;
                            ok=4;
                        }
                    }
                    else if (ok==4)
                    {
                        udep=udepnou;
                        vdep=vdepnou;
                        uv[0]=udep;
                        uv[1]=vdep;
                        if ((udep+1e-10>umax) && (vdep+1e-10>vmax))
                            ok=-1;
                        else
                            ok=0;
                    }
                }
                while (ok<5);
            }
 
        }
    }
    //insertion des points de controles
    for (int i=0;i<nb_pt_ctrl;i++)
    {
        double xyz[3];
        xyz[0]=point_controle[4*i];
        xyz[1]=point_controle[4*i+1];
        xyz[2]=point_controle[4*i+2];
        double dn=point_controle[4*i+3];
        double dist=zonezoom*(dng-dn);
        BOITE_3D boite(xyz[0]-dist,xyz[1]-dist,xyz[2]-dist,xyz[0]+dist,xyz[1]+dist,xyz[2]+dist);
        DOUBLEN<1> ecart(&dn);
        CT_POINT<1>* pt=new CT_POINT<1>(xyz[0],xyz[1],xyz[2],ecart,boite);
        lstpoint.ajouter(pt);
    }
}
 
 
 
 
void FCT_GENERATEUR_FRONTIERE::construit(int pas_echan,int pas_estima,MG_GROUPE_TOPOLOGIQUE* mggt)
{
    pas_echantillon=pas_echan;
    pas_estimation=pas_estima;
    echantillonner(mggt);
    estimation();
}
 
 
#pragma package(smart_init)