Back to index

salome-geom  6.5.0
GEOM_WireframeFace.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2007-2012  CEA/DEN, EDF R&D, OPEN CASCADE
00002 //
00003 // This library is free software; you can redistribute it and/or
00004 // modify it under the terms of the GNU Lesser General Public
00005 // License as published by the Free Software Foundation; either
00006 // version 2.1 of the License.
00007 //
00008 // This library is distributed in the hope that it will be useful,
00009 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00010 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00011 // Lesser General Public License for more details.
00012 //
00013 // You should have received a copy of the GNU Lesser General Public
00014 // License along with this library; if not, write to the Free Software
00015 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00016 //
00017 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00018 //
00019 
00020 #include "GEOM_WireframeFace.h" 
00021  
00022 #include <vtkObjectFactory.h> 
00023  
00024 #include <vtkPoints.h> 
00025 #include <vtkCellArray.h> 
00026 
00027 #include <vtkPolyDataMapper.h>  
00028 #include <vtkPolyData.h>  
00029  
00030 #include <Precision.hxx>
00031 #include <BRepTools.hxx>
00032 #include <TopExp_Explorer.hxx>
00033 #include <Geom2dHatch_Hatcher.hxx>
00034 #include <Geom2dHatch_Intersector.hxx>
00035 #include <TColStd_Array1OfReal.hxx>
00036 #include <TColStd_Array1OfInteger.hxx>
00037  
00038 #include <TopoDS.hxx> 
00039 #include <TopoDS_Edge.hxx> 
00040 #include <BRep_Tool.hxx>
00041 #include <Geom2d_TrimmedCurve.hxx>
00042 #include <Geom2d_Line.hxx>
00043 #include <gp_Dir2d.hxx>
00044 #include <gp_Pnt2d.hxx>
00045  
00046 #include <Geom2dHatch_Hatcher.hxx>
00047 #include <HatchGen_Domain.hxx>
00048 
00049 #include <Adaptor3d_HCurve.hxx>
00050 
00051 vtkStandardNewMacro(GEOM_WireframeFace);
00052  
00053 GEOM_WireframeFace::GEOM_WireframeFace(): 
00054   Discret(15)
00055 { 
00056   NbIso[0] = 1;
00057   NbIso[1] = 1;
00058 } 
00059  
00060 GEOM_WireframeFace::~GEOM_WireframeFace() 
00061 { 
00062 } 
00063  
00064 void
00065 GEOM_WireframeFace:: 
00066 Execute()
00067 {
00068   vtkPolyData* aPolyData = GetOutput();
00069   aPolyData->Allocate();
00070   vtkPoints* aPts = vtkPoints::New();
00071   aPolyData->SetPoints(aPts);
00072   aPts->Delete();
00073 
00074   TFaceSet::Iterator anIter(myFaceSet);
00075   for(; anIter.More(); anIter.Next()){
00076     const TopoDS_Face& aFace = anIter.Value();
00077     OCC2VTK(aFace,aPolyData,aPts,NbIso,Discret);
00078   }
00079 }
00080 
00081 void GEOM_WireframeFace::SetNbIso(const int theNb[2])
00082 {
00083   if ( theNb[0] == NbIso[0] && theNb[1] == NbIso[1])
00084     return;
00085 
00086   NbIso[0] = theNb[0];
00087   NbIso[1] = theNb[1];
00088 
00089   Modified();
00090 }
00091 
00092 void GEOM_WireframeFace::GetNbIso(int &theNbU,int &theNbV)
00093 {
00094   theNbU = NbIso[0];
00095   theNbV = NbIso[1];
00096 }
00097 
00098 void  
00099 GEOM_WireframeFace:: 
00100 OCC2VTK(const TopoDS_Face& theFace,
00101         vtkPolyData* thePolyData,
00102                     vtkPoints* thePts,  
00103         const int theNbIso[2], 
00104         const int theDiscret) 
00105 { 
00106   TopoDS_Face aFace = theFace; 
00107   aFace.Orientation(TopAbs_FORWARD);
00108   CreateIso(aFace,theNbIso,theDiscret,thePolyData,thePts);
00109 }
00110 
00111 void 
00112 GEOM_WireframeFace:: 
00113 CreateIso(const TopoDS_Face& theFace,
00114           const int theNbIso[2], 
00115           const int theDiscret, 
00116           vtkPolyData* thePolyData,
00117           vtkPoints* thePts)
00118 {
00119   // Constants for iso building
00120   static Standard_Real INTERSECTOR_CONFUSION = 1.e-10 ; // -8 ;
00121   static Standard_Real INTERSECTOR_TANGENCY  = 1.e-10 ; // -8 ;
00122 
00123   static Standard_Real HATHCER_CONFUSION_2D = 1.e-8 ;
00124   static Standard_Real HATHCER_CONFUSION_3D = 1.e-8 ;
00125 
00126   Geom2dHatch_Hatcher 
00127     aHatcher(Geom2dHatch_Intersector(INTERSECTOR_CONFUSION,
00128                                      INTERSECTOR_TANGENCY),
00129                          HATHCER_CONFUSION_2D,
00130                          HATHCER_CONFUSION_3D,
00131                                      Standard_True,
00132                                      Standard_False);
00133   
00134   Standard_Real anUMin, anUMax, aVMin, aVMax;
00135   TColStd_Array1OfReal anUPrm(0, theNbIso[0]), aVPrm(0, theNbIso[1]);
00136   TColStd_Array1OfInteger anUInd(0, theNbIso[0]), aVInd(0, theNbIso[1]);
00137 
00138   anUInd.Init(0);
00139   aVInd.Init(0);
00140 
00141   //-----------------------------------------------------------------------
00142   // If the Min Max bounds are infinite, there are bounded to Infinite
00143   // value.
00144   //-----------------------------------------------------------------------
00145   BRepTools::UVBounds(theFace, anUMin, anUMax, aVMin, aVMax) ;
00146   Standard_Boolean InfiniteUMin = Precision::IsNegativeInfinite (anUMin) ;
00147   Standard_Boolean InfiniteUMax = Precision::IsPositiveInfinite (anUMax) ;
00148   Standard_Boolean InfiniteVMin = Precision::IsNegativeInfinite (aVMin) ;
00149   Standard_Boolean InfiniteVMax = Precision::IsPositiveInfinite (aVMax) ;
00150 
00151   static float VTKINFINITE = 1.0E38;
00152   if(InfiniteUMin && InfiniteUMax){
00153     anUMin = - VTKINFINITE ;
00154     anUMax =   VTKINFINITE ;
00155   }else if(InfiniteUMin){
00156     anUMin = anUMax - VTKINFINITE ;
00157   }else if(InfiniteUMax){
00158     anUMax = anUMin + VTKINFINITE ;
00159   }
00160 
00161   if(InfiniteVMin && InfiniteVMax){
00162     aVMin = - VTKINFINITE ;
00163     aVMax =   VTKINFINITE ;
00164   }else if(InfiniteVMin){
00165     aVMin = aVMax - VTKINFINITE ;
00166   }else if(InfiniteVMax){
00167     aVMax = aVMin + VTKINFINITE ;
00168   }
00169 
00170   //-----------------------------------------------------------------------
00171   // Retreiving the edges and loading them into the hatcher.
00172   //-----------------------------------------------------------------------
00173   TopExp_Explorer ExpEdges(theFace, TopAbs_EDGE);
00174   for(; ExpEdges.More(); ExpEdges.Next()){
00175     const TopoDS_Edge& anEdge = TopoDS::Edge(ExpEdges.Current());
00176     Standard_Real U1, U2 ;
00177     const Handle(Geom2d_Curve) PCurve = 
00178       BRep_Tool::CurveOnSurface(anEdge, theFace, U1, U2) ;
00179 
00180     if(PCurve.IsNull() || U1 == U2)
00181       return;
00182 
00183     //-- Test if a TrimmedCurve is necessary
00184     if(Abs(PCurve->FirstParameter()-U1) <= Precision::PConfusion() &&
00185              Abs(PCurve->LastParameter()-U2) <= Precision::PConfusion())
00186     { 
00187       aHatcher.AddElement(PCurve, anEdge.Orientation()) ;      
00188     }else{ 
00189       if(!PCurve->IsPeriodic()){
00190               Handle(Geom2d_TrimmedCurve) TrimPCurve =
00191           Handle(Geom2d_TrimmedCurve)::DownCast(PCurve);
00192               if(!TrimPCurve.IsNull()){
00193           Handle_Geom2d_Curve aBasisCurve = TrimPCurve->BasisCurve();
00194                 if(aBasisCurve->FirstParameter()-U1 > Precision::PConfusion() ||
00195                    U2-aBasisCurve->LastParameter() > Precision::PConfusion()) 
00196           {
00197                   aHatcher.AddElement(PCurve, anEdge.Orientation()) ;      
00198                   return;
00199                 }
00200               }else{
00201                 if(PCurve->FirstParameter()-U1 > Precision::PConfusion()){
00202                   U1=PCurve->FirstParameter();
00203                 }
00204                 if(U2-PCurve->LastParameter()  > Precision::PConfusion()){
00205                   U2=PCurve->LastParameter();
00206                 }
00207               }
00208       }
00209       Handle(Geom2d_TrimmedCurve) TrimPCurve = 
00210         new Geom2d_TrimmedCurve(PCurve, U1, U2);
00211       aHatcher.AddElement(TrimPCurve, anEdge.Orientation());
00212     }
00213   }
00214 
00215 
00216   //-----------------------------------------------------------------------
00217   // Loading and trimming the hatchings.
00218   //-----------------------------------------------------------------------
00219   Standard_Integer IIso;
00220   Standard_Real DeltaU = Abs(anUMax - anUMin) ;
00221   Standard_Real DeltaV = Abs(aVMax - aVMin) ;
00222   Standard_Real confusion = Min(DeltaU, DeltaV) * HATHCER_CONFUSION_3D ;
00223   aHatcher.Confusion3d (confusion) ;
00224 
00225   if ( theNbIso[0] ) {
00226     Standard_Real StepU = DeltaU / (Standard_Real)theNbIso[0];
00227     if(StepU > confusion){
00228       Standard_Real UPrm = anUMin + StepU / 2.;
00229       gp_Dir2d Dir(0., 1.) ;
00230       for(IIso = 1 ; IIso <= theNbIso[0] ; IIso++) {
00231         anUPrm(IIso) = UPrm ;
00232         gp_Pnt2d Ori (UPrm, 0.) ;
00233         Geom2dAdaptor_Curve HCur (new Geom2d_Line (Ori, Dir)) ;
00234         anUInd(IIso) = aHatcher.AddHatching (HCur) ;
00235         UPrm += StepU ;
00236       }
00237     }
00238   }
00239 
00240   if ( theNbIso[1] ) {
00241     Standard_Real StepV = DeltaV / (Standard_Real) theNbIso[1] ;
00242     if(StepV > confusion){
00243       Standard_Real VPrm = aVMin + StepV / 2.;
00244       gp_Dir2d Dir(1., 0.);
00245       for(IIso = 1 ; IIso <= theNbIso[1] ; IIso++){
00246         aVPrm(IIso) = VPrm;
00247         gp_Pnt2d Ori (0., VPrm);
00248         Geom2dAdaptor_Curve HCur(new Geom2d_Line (Ori, Dir));
00249         aVInd(IIso) = aHatcher.AddHatching (HCur) ;
00250         VPrm += StepV ;
00251       }
00252     }
00253   }
00254 
00255   //-----------------------------------------------------------------------
00256   // Computation.
00257   //-----------------------------------------------------------------------
00258   aHatcher.Trim() ;
00259 
00260   Standard_Integer aNbDom = 0 ; // for debug purpose
00261   Standard_Integer Index ;
00262 
00263   for(IIso = 1 ; IIso <= theNbIso[0] ; IIso++){
00264     Index = anUInd(IIso) ;
00265     if(Index != 0){
00266       if(aHatcher.TrimDone(Index) && !aHatcher.TrimFailed(Index)){
00267               aHatcher.ComputeDomains(Index);
00268               if(aHatcher.IsDone (Index)) 
00269           aNbDom = aHatcher.NbDomains (Index);
00270       }
00271     }
00272   }
00273 
00274   for(IIso = 1 ; IIso <= theNbIso[1] ; IIso++){
00275     Index = aVInd(IIso);
00276     if(Index != 0){
00277       if(aHatcher.TrimDone (Index) && !aHatcher.TrimFailed(Index)){
00278               aHatcher.ComputeDomains (Index);
00279               if(aHatcher.IsDone (Index)) 
00280           aNbDom = aHatcher.NbDomains (Index);
00281       }
00282     }
00283   }
00284 
00285   //-----------------------------------------------------------------------
00286   // Push iso lines in vtk kernel
00287   //-----------------------------------------------------------------------
00288   for(Standard_Integer UIso = anUPrm.Lower() ; UIso <= anUPrm.Upper(); UIso++){
00289     Standard_Integer UInd = anUInd.Value(UIso);
00290     if(UInd != 0){
00291       Standard_Real UPrm = anUPrm.Value(UIso);
00292       if(aHatcher.IsDone(UInd)){
00293               Standard_Integer NbDom = aHatcher.NbDomains(UInd);
00294               for(Standard_Integer IDom = 1 ; IDom <= NbDom ; IDom++){
00295                 const HatchGen_Domain& Dom = aHatcher.Domain(UInd, IDom) ;
00296                 Standard_Real V1 = Dom.HasFirstPoint()? Dom.FirstPoint().Parameter(): aVMin - VTKINFINITE;
00297                 Standard_Real V2 = Dom.HasSecondPoint()? Dom.SecondPoint().Parameter(): aVMax + VTKINFINITE;
00298                 CreateIso_(theFace, GeomAbs_IsoU, UPrm, V1, V2, theDiscret, thePolyData, thePts);
00299         }
00300             }
00301     }
00302   }
00303 
00304   for(Standard_Integer VIso = aVPrm.Lower() ; VIso <= aVPrm.Upper(); VIso++){
00305     Standard_Integer VInd = aVInd.Value(VIso);
00306     if(VInd != 0){
00307       Standard_Real VPrm = aVPrm.Value(VIso);
00308       if(aHatcher.IsDone (VInd)){
00309               Standard_Integer NbDom = aHatcher.NbDomains(VInd);
00310               for (Standard_Integer IDom = 1 ; IDom <= NbDom ; IDom++){
00311                 const HatchGen_Domain& Dom = aHatcher.Domain(VInd, IDom);
00312                 Standard_Real U1 = Dom.HasFirstPoint()? Dom.FirstPoint().Parameter(): aVMin - VTKINFINITE;
00313                 Standard_Real U2 = Dom.HasSecondPoint()? Dom.SecondPoint().Parameter(): aVMax + VTKINFINITE;
00314             CreateIso_(theFace, GeomAbs_IsoV, VPrm, U1, U2, theDiscret, thePolyData, thePts);
00315               }
00316       }
00317     }
00318   }
00319 }
00320 
00321  
00322 
00323 void 
00324 GEOM_WireframeFace:: 
00325 CreateIso_(const TopoDS_Face& theFace,
00326            GeomAbs_IsoType theIsoType, 
00327            Standard_Real Par, 
00328            Standard_Real T1,
00329            Standard_Real T2,
00330            const int theDiscret, 
00331            vtkPolyData* thePolyData,
00332            vtkPoints* thePts)
00333 {
00334   Standard_Real U1, U2, V1, V2, stepU=0., stepV=0.;
00335   Standard_Integer j;
00336   gp_Pnt P;
00337 
00338   TopLoc_Location aLoc;
00339   const Handle(Geom_Surface)& S = BRep_Tool::Surface(theFace,aLoc);
00340 
00341   if(!S.IsNull()){
00342     BRepAdaptor_Surface S(theFace,Standard_False);
00343       
00344     GeomAbs_SurfaceType SurfType = S.GetType();
00345 
00346     GeomAbs_CurveType CurvType = GeomAbs_OtherCurve;
00347 
00348     Standard_Integer Intrv, nbIntv;
00349     Standard_Integer nbUIntv = S.NbUIntervals(GeomAbs_CN);
00350     Standard_Integer nbVIntv = S.NbVIntervals(GeomAbs_CN);
00351     TColStd_Array1OfReal TI(1,Max(nbUIntv, nbVIntv)+1);
00352 
00353     if(theIsoType == GeomAbs_IsoU){
00354       S.VIntervals(TI, GeomAbs_CN);
00355       V1 = Max(T1, TI(1));
00356       V2 = Min(T2, TI(2));
00357       U1 = Par;
00358       U2 = Par;
00359       stepU = 0;
00360       nbIntv = nbVIntv;
00361     }else{
00362       S.UIntervals(TI, GeomAbs_CN);
00363       U1 = Max(T1, TI(1));
00364       U2 = Min(T2, TI(2));
00365       V1 = Par;
00366       V2 = Par;
00367       stepV = 0;
00368       nbIntv = nbUIntv;
00369     }   
00370         
00371     S.D0(U1,V1,P);
00372     MoveTo(P,thePts);
00373 
00374     for(Intrv = 1; Intrv <= nbIntv; Intrv++){
00375       if(TI(Intrv) <= T1 && TI(Intrv + 1) <= T1)
00376         continue;
00377       if(TI(Intrv) >= T2 && TI(Intrv + 1) >= T2)
00378               continue;
00379       if(theIsoType == GeomAbs_IsoU){
00380               V1 = Max(T1, TI(Intrv));
00381               V2 = Min(T2, TI(Intrv + 1));
00382               stepV = (V2 - V1) / theDiscret;
00383       }else{
00384               U1 = Max(T1, TI(Intrv));
00385               U2 = Min(T2, TI(Intrv + 1));
00386               stepU = (U2 - U1) / theDiscret;
00387       }
00388 
00389       switch (SurfType) {
00390       case GeomAbs_Plane :
00391               break;
00392       case GeomAbs_Cylinder :
00393       case GeomAbs_Cone :
00394         if(theIsoType == GeomAbs_IsoV){
00395                 for(j = 1; j < theDiscret; j++){
00396                   U1 += stepU;
00397                   V1 += stepV;
00398                   S.D0(U1,V1,P);
00399                   DrawTo(P,thePolyData,thePts);
00400                 }
00401               }
00402               break;
00403       case GeomAbs_Sphere :
00404       case GeomAbs_Torus :
00405       case GeomAbs_OffsetSurface :
00406       case GeomAbs_OtherSurface :
00407         for(j = 1; j < theDiscret; j++){
00408                 U1 += stepU;
00409                 V1 += stepV;
00410                 S.D0(U1,V1,P);
00411                 DrawTo(P,thePolyData,thePts);
00412               }
00413               break;
00414       case GeomAbs_BezierSurface :
00415       case GeomAbs_BSplineSurface :
00416         for(j = 1; j <= theDiscret/2; j++){
00417           Standard_Real aStep = (theIsoType == GeomAbs_IsoV) ? stepU*2. : stepV*2.;
00418                 CreateIso__(S, theIsoType, U1, V1, aStep, thePolyData, thePts);
00419                 U1 += stepU*2.;
00420                 V1 += stepV*2.;
00421               }
00422               break;
00423       case GeomAbs_SurfaceOfExtrusion :
00424       case GeomAbs_SurfaceOfRevolution :
00425         if((theIsoType == GeomAbs_IsoV && SurfType == GeomAbs_SurfaceOfRevolution) ||
00426                  (theIsoType == GeomAbs_IsoU && SurfType == GeomAbs_SurfaceOfExtrusion)) 
00427         {
00428                 if(SurfType == GeomAbs_SurfaceOfExtrusion) 
00429             break;
00430                 for(j = 1; j < theDiscret; j++){
00431                   U1 += stepU;
00432                   V1 += stepV;
00433                   S.D0(U1,V1,P);
00434                   DrawTo(P,thePolyData,thePts);
00435                 }
00436               }else{
00437                 CurvType = (S.BasisCurve())->GetType();
00438                 switch(CurvType){
00439                 case GeomAbs_Line :
00440                   break;
00441                 case GeomAbs_Circle :
00442                 case GeomAbs_Ellipse :
00443                   for (j = 1; j < theDiscret; j++) {
00444                     U1 += stepU;
00445                     V1 += stepV;
00446                     S.D0(U1,V1,P);
00447                     DrawTo(P,thePolyData,thePts);
00448                   }
00449                   break;
00450                 case GeomAbs_Parabola :
00451                 case GeomAbs_Hyperbola :
00452                 case GeomAbs_BezierCurve :
00453                 case GeomAbs_BSplineCurve :
00454                 case GeomAbs_OtherCurve :
00455                   for(j = 1; j <= theDiscret/2; j++){
00456               Standard_Real aStep = (theIsoType == GeomAbs_IsoV) ? stepU*2. : stepV*2.;
00457                   CreateIso__(S, theIsoType, U1, V1, aStep, thePolyData, thePts);
00458                     U1 += stepU*2.;
00459                     V1 += stepV*2.;
00460                   }
00461                   break;
00462                 }
00463               }
00464       }
00465     }
00466     S.D0(U2,V2,P);
00467     DrawTo(P,thePolyData,thePts);
00468   }  
00469 }
00470  
00471  
00472  
00473  
00474 void 
00475 GEOM_WireframeFace:: 
00476 CreateIso__(const BRepAdaptor_Surface& theSurface, 
00477             GeomAbs_IsoType theIsoType,
00478                                     Standard_Real& theU, 
00479                                     Standard_Real& theV, 
00480                                     Standard_Real theStep, 
00481             vtkPolyData* thePolyData,
00482             vtkPoints* thePts)
00483 {
00484   gp_Pnt Pl, Pr, Pm;
00485   if (theIsoType == GeomAbs_IsoU) {
00486     theSurface.D0(theU, theV, Pl);
00487     theSurface.D0(theU, theV + theStep/2., Pm);
00488     theSurface.D0(theU, theV + theStep, Pr);
00489   } else {
00490     theSurface.D0(theU, theV, Pl);
00491     theSurface.D0(theU + theStep/2., theV, Pm);
00492     theSurface.D0(theU + theStep, theV, Pr);
00493   }
00494 
00495   static Standard_Real ISO_RATIO = 1.001;
00496   if (Pm.Distance(Pl) + Pm.Distance(Pr) <= ISO_RATIO*Pl.Distance(Pr)) {
00497     DrawTo(Pr,thePolyData,thePts);
00498   } else {
00499     if (theIsoType == GeomAbs_IsoU) {
00500       CreateIso__(theSurface, theIsoType, theU, theV, theStep/2, thePolyData, thePts);
00501       Standard_Real aLocalV = theV + theStep/2 ;
00502       CreateIso__(theSurface, theIsoType, theU, aLocalV , theStep/2, thePolyData, thePts);
00503     } else {
00504       CreateIso__(theSurface, theIsoType, theU, theV, theStep/2, thePolyData, thePts);
00505       Standard_Real aLocalU = theU + theStep/2 ;
00506       CreateIso__(theSurface, theIsoType, aLocalU , theV, theStep/2, thePolyData, thePts);
00507     }
00508   }
00509 }