Back to index

salome-gui  6.5.0
vtkEDFCutter.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2010-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 "vtkEDFCutter.h"
00021 
00022 #include "vtkInformationVector.h"
00023 #include "vtkInformation.h"
00024 #include "vtkSmartPointer.h"
00025 #include "vtkGenericCell.h"
00026 #include "vtkPolyData.h"
00027 #include "vtkObjectFactory.h"
00028 #include "vtkIdTypeArray.h"
00029 #include "vtkCellData.h"
00030 #include "vtkCellArray.h"
00031 #include "vtkIdList.h"
00032 
00033 #include <list>
00034 #include <set>
00035 #include <map>
00036 #include <deque>
00037 
00038 class vtkEDFEdge
00039 {
00040 public :
00041   vtkIdType v0;
00042   vtkIdType v1;
00043 
00044   vtkEDFEdge(vtkIdType a, vtkIdType b) : v0(a), v1(b){}
00045   vtkEDFEdge(){}
00046 };
00047 
00048 bool operator == (const vtkEDFEdge& e0, const vtkEDFEdge& e1)
00049 {
00050   return (e0.v0 == e1.v0 && e0.v1 == e1.v1) ||
00051       (e0.v1 == e1.v0 && e0.v0 == e1.v1);
00052 }
00053 
00054 bool operator != (const vtkEDFEdge& e0, const vtkEDFEdge& e1)
00055 {
00056   return !(e0==e1);
00057 }
00058 
00059 bool operator < (const vtkEDFEdge& e0, const vtkEDFEdge& e1)
00060 {
00061   vtkEDFEdge the_e0;
00062   vtkEDFEdge the_e1;
00063   if(e0.v0 < e0.v1)
00064     {
00065     the_e0.v0 = e0.v0;
00066     the_e0.v1 = e0.v1;
00067     }
00068   else
00069     {
00070     the_e0.v0 = e0.v1;
00071     the_e0.v1 = e0.v0;
00072     }
00073   if(e1.v0 < e1.v1)
00074     {
00075     the_e1.v0 = e1.v0;
00076     the_e1.v1 = e1.v1;
00077     }
00078   else
00079     {
00080     the_e1.v0 = e1.v1;
00081     the_e1.v1 = e1.v0;
00082     }
00083 
00084   if(the_e0.v0 == the_e1.v0)
00085     return (the_e0.v1 < the_e1.v1);
00086 
00087   return the_e0.v0 < the_e1.v0;
00088 }
00089 
00090 vtkStandardNewMacro(vtkEDFCutter);
00091 vtkCxxRevisionMacro(vtkEDFCutter, "0.0");
00092 
00093 vtkEDFCutter::vtkEDFCutter()
00094 {
00095   this->OriginalCellIdsName = NULL;
00096 }
00097 
00098 vtkEDFCutter::~vtkEDFCutter()
00099 {
00100   this->SetOriginalCellIdsName(NULL);
00101 }
00102 
00103 int vtkEDFCutter::RequestData(vtkInformation * request,
00104                         vtkInformationVector ** inVector,
00105                         vtkInformationVector * outVector)
00106 {
00107   // get the info objects
00108   vtkInformation *inInfo = inVector[0]->GetInformationObject(0);
00109   vtkInformation *outInfo = outVector->GetInformationObject(0);
00110 
00111   // get the input and output
00112   vtkDataSet *input = vtkDataSet::SafeDownCast(
00113     inInfo->Get(vtkDataObject::DATA_OBJECT()));
00114 
00115   vtkSmartPointer<vtkIdTypeArray> cellIdArray =
00116       vtkSmartPointer<vtkIdTypeArray>::New();
00117   cellIdArray->SetName(this->GetOriginalCellIdsName());
00118   cellIdArray->SetNumberOfComponents(1);
00119   cellIdArray->SetNumberOfTuples(input->GetNumberOfCells());
00120   for(vtkIdType id=0; id < cellIdArray->GetNumberOfTuples(); id++)
00121     {
00122     cellIdArray->SetTuple1(id, id);
00123     }
00124   input->GetCellData()->AddArray(cellIdArray);
00125 
00126   int ret = this->Superclass::RequestData(request, inVector, outVector);
00127 
00128   if(ret == 0)
00129     return 0;
00130 
00131   vtkPolyData *output = vtkPolyData::SafeDownCast(
00132     outInfo->Get(vtkDataObject::DATA_OBJECT()));
00133 
00134   vtkSmartPointer<vtkPolyData> tmpOutput;
00135   tmpOutput.TakeReference(output->NewInstance());
00136 
00137   this->AssembleOutputTriangles(output, tmpOutput);
00138 
00139   output->ShallowCopy(tmpOutput);
00140 
00141   return ret;
00142 }
00143 
00144 
00145 void  vtkEDFCutter::AssembleOutputTriangles(vtkPolyData* inpd,
00146                                             vtkPolyData* outpd)
00147 {
00148   outpd->ShallowCopy(inpd);
00149 
00150   vtkIdTypeArray* originalCellIds = vtkIdTypeArray::SafeDownCast(
00151       inpd->GetCellData()->GetArray(this->GetOriginalCellIdsName()));
00152 
00153   if(originalCellIds == NULL)
00154     {
00155     return;
00156     }
00157 
00158   outpd->GetCellData()->Initialize();
00159   outpd->GetCellData()->CopyAllocate(inpd->GetCellData());
00160 
00161   vtkSmartPointer<vtkCellArray> verts = vtkSmartPointer<vtkCellArray>::New();
00162   vtkSmartPointer<vtkCellArray> lines = vtkSmartPointer<vtkCellArray>::New();
00163   vtkSmartPointer<vtkCellArray> polys = vtkSmartPointer<vtkCellArray>::New();
00164   vtkSmartPointer<vtkCellArray> strips = vtkSmartPointer<vtkCellArray>::New();
00165   outpd->SetVerts(verts);
00166   outpd->SetLines(lines);
00167   outpd->SetPolys(polys);
00168   outpd->SetStrips(strips);
00169 
00170   for(vtkIdType cellId=0; cellId<inpd->GetNumberOfCells(); cellId++)
00171     {
00172     unsigned char type = inpd->GetCellType(cellId);
00173     if(type != VTK_TRIANGLE)
00174       {
00175       vtkIdType npts;
00176       vtkIdType* pts = NULL;
00177       inpd->GetCellPoints(cellId, npts, pts);
00178       vtkIdType newCellId =
00179           outpd->InsertNextCell(type, npts, pts);
00180       outpd->GetCellData()->CopyData(inpd->GetCellData(), cellId, newCellId);
00181       }
00182     else
00183       {
00184       vtkIdType cellIdEnd = cellId+1;
00185       vtkIdType originalCellId = originalCellIds->GetValue(cellId);
00186       while(cellIdEnd < inpd->GetNumberOfCells() &&
00187             inpd->GetCellType(cellIdEnd) == VTK_TRIANGLE &&
00188             originalCellIds->GetValue(cellIdEnd) == originalCellId)
00189         {
00190         cellIdEnd++;
00191         }
00192 
00193       // all triangles from cellId to cellIdEnd come from the same
00194       // original cell.
00195 
00196       // A batch is composed of triangles which are connected by the edges.
00197       std::map<vtkIdType, std::set<vtkIdType> > connectedTriangles;
00198       for(vtkIdType firstCell = cellId; firstCell < cellIdEnd-1; firstCell++)
00199         {
00200         vtkIdType npts;
00201         vtkIdType* pts = NULL;
00202         inpd->GetCellPoints(firstCell, npts, pts);
00203         vtkEDFEdge fe0 = vtkEDFEdge(pts[0], pts[1]);
00204         vtkEDFEdge fe1 = vtkEDFEdge(pts[1], pts[2]);
00205         vtkEDFEdge fe2 = vtkEDFEdge(pts[2], pts[0]);
00206         for(vtkIdType secondCell = firstCell+1; secondCell < cellIdEnd; secondCell++)
00207           {
00208           vtkIdType snpts;
00209           vtkIdType* spts = NULL;
00210           inpd->GetCellPoints(secondCell, snpts, spts);
00211           vtkEDFEdge se0 = vtkEDFEdge(spts[0], spts[1]);
00212           vtkEDFEdge se1 = vtkEDFEdge(spts[1], spts[2]);
00213           vtkEDFEdge se2 = vtkEDFEdge(spts[2], spts[0]);
00214 
00215           if(fe0 == se0 || fe0 == se1 || fe0 == se2 ||
00216              fe1 == se0 || fe1 == se1 || fe1 == se2 ||
00217              fe2 == se0 || fe2 == se1 || fe2 == se2)
00218             {
00219             connectedTriangles[firstCell].insert(secondCell);
00220             connectedTriangles[secondCell].insert(firstCell);
00221             }
00222           }
00223         }
00224 
00225       std::set<vtkIdType> visitedCell;
00226       for(vtkIdType id=cellId; id<cellIdEnd; id++)
00227         {
00228         if(visitedCell.find(id) != visitedCell.end())
00229           continue;
00230 
00231         // if this cell has not yet been visited, I create a batch of all
00232         // cells connected to this one
00233 
00234         visitedCell.insert(id);
00235         std::set<vtkIdType> batch;
00236         std::list<vtkIdType> cellList;
00237         cellList.push_back(id);
00238         while(cellList.size() > 0)
00239           {
00240           vtkIdType currentId = *(cellList.begin());
00241           batch.insert(currentId);
00242           cellList.pop_front();
00243           if(connectedTriangles.find(currentId) != connectedTriangles.end())
00244             {
00245             const std::set<vtkIdType>& adj = connectedTriangles[currentId];
00246             std::set<vtkIdType>::const_iterator it = adj.begin();
00247             while(it != adj.end())
00248               {
00249               vtkIdType other = *it;
00250               if(visitedCell.find(other) == visitedCell.end())
00251                 {
00252                 cellList.push_back(other);
00253                 visitedCell.insert(other);
00254                 }
00255               it++;
00256               }
00257             }
00258           }
00259 
00260 
00261 
00262         // then I add this batch to the output,
00263         // creating a unique cell for the whole batch.
00264 
00265         if(batch.size() == 1)
00266           {
00267           vtkIdType tid = *(batch.begin());
00268           vtkIdType npts;
00269           vtkIdType* pts = NULL;
00270           inpd->GetCellPoints(tid, npts, pts);
00271           vtkIdType newCellId =
00272               outpd->InsertNextCell(VTK_TRIANGLE, npts, pts);
00273           outpd->GetCellData()->CopyData(inpd->GetCellData(), cellId, newCellId);
00274           }
00275         else if(batch.size() == 2)
00276           { // two triangles connected together --> create a VTK_QUAD
00277           vtkIdType fid = *(batch.begin());
00278           vtkIdType sid = *(batch.rbegin());
00279           vtkIdType fnpts;
00280           vtkIdType* fpts = NULL;
00281           inpd->GetCellPoints(fid, fnpts, fpts);
00282           vtkIdType snpts;
00283           vtkIdType* spts = NULL;
00284           inpd->GetCellPoints(sid, snpts, spts);
00285 
00286           int findex = 0;
00287           vtkIdType fv = fpts[findex];
00288           while(((fv == spts[0]) ||
00289                  (fv == spts[1]) ||
00290                  (fv == spts[2])) && findex < 3)
00291             {
00292             findex++;
00293             fv = fpts[findex];
00294             }
00295           if(findex == 3)
00296             {// this is a degenerate case : one of the triangles use
00297             // only 2 vertices
00298             findex = 0;
00299             }
00300           int sindex = 0;
00301           vtkIdType sv = spts[sindex];
00302           while(((sv == fpts[0]) ||
00303                  (sv == fpts[1]) ||
00304                  (sv == fpts[2])) && sindex < 3)
00305             {
00306             sindex++;
00307             sv = spts[sindex];
00308             }
00309           if(sindex == 3)
00310             {// this is a degenerate case : one of the triangles use
00311             // only 2 vertices
00312             sindex = 0;
00313             }
00314 
00315           vtkIdType pts[4];
00316           pts[0] = fpts[findex];
00317           pts[1] = fpts[(findex+1)%3];
00318           pts[2] = spts[sindex];
00319           pts[3] = fpts[(findex+2)%3];
00320 
00321           vtkIdType newCellId = outpd->InsertNextCell(VTK_QUAD, 4, pts);
00322           outpd->GetCellData()->CopyData(inpd->GetCellData(), cellId, newCellId);
00323           }
00324         else
00325           {
00326           std::deque<vtkEDFEdge> contour;
00327 
00328           std::list<vtkIdType> toVisit;
00329           std::set<vtkIdType> visited;
00330 
00331           toVisit.push_back(*(batch.begin()));
00332 
00333           std::set<vtkIdType> triedAgain;
00334 
00335           while(toVisit.size() > 0)
00336             {
00337             vtkIdType currentId = *(toVisit.begin());
00338             toVisit.pop_front();
00339             if(visited.find(currentId) != visited.end())
00340               continue;
00341 
00342             visited.insert(currentId);
00343             const std::set<vtkIdType>& adj = connectedTriangles[currentId];
00344             std::set<vtkIdType>::const_iterator it = adj.begin();
00345             while(it != adj.end())
00346               {
00347               vtkIdType adjid = *it;
00348               it++;
00349               if(visited.find(adjid) != visited.end())
00350                 continue;
00351 
00352               toVisit.push_back(adjid);
00353               }
00354 
00355             vtkIdType npts;
00356             vtkIdType* pts = NULL;
00357             inpd->GetCellPoints(currentId, npts, pts);
00358             vtkEDFEdge edge[3] = {vtkEDFEdge(pts[0], pts[1]),
00359                                   vtkEDFEdge(pts[1], pts[2]),
00360                                   vtkEDFEdge(pts[2], pts[0])};
00361 
00362             // special case : initialization of the contour
00363             if(contour.size() == 0)
00364               {
00365               contour.push_back(edge[0]);
00366               contour.push_back(edge[1]);
00367               contour.push_back(edge[2]);
00368               continue;
00369               }
00370 
00371             // Find which edge of the contour
00372             // is connected to the current triangle
00373             int orient = 0;
00374             std::deque<vtkEDFEdge>::iterator contourit = contour.begin();
00375             bool found = false;
00376             while(contourit != contour.end())
00377               {
00378               vtkEDFEdge& e = *contourit;
00379               for(orient = 0; orient<3; orient++)
00380                 {
00381                 if(e == edge[orient])
00382                   {
00383                   found = true;
00384                   break;
00385                   }
00386                 }
00387               if(found)
00388                 break;
00389 
00390               contourit++;
00391               }
00392             if(contourit == contour.end())
00393               {// this triangle is not connected to the current contour.
00394               // put it back in the queue for later processing
00395               if(triedAgain.find(currentId) == triedAgain.end())
00396                 {
00397                 triedAgain.insert(currentId);
00398                 toVisit.push_back(currentId);
00399                 visited.erase(currentId);
00400                 continue;
00401                 }
00402               else
00403                 {
00404                 vtkDebugMacro( << "triangle " << currentId
00405                   << "could not be added to the contour of the current batch");
00406                 continue;
00407                 }
00408               }
00409             // if I reach this point, I will add the triangle to the contour
00410             // --> the contour will be modified and I can try again
00411             // to add the previously rejected triangles
00412             triedAgain.clear();
00413 
00414             // Now, merge the edges of the current triangle with
00415             // the contour
00416             vtkEDFEdge& tbeforeedge = edge[(orient+1)%3];
00417             vtkEDFEdge& tafteredge = edge[(orient+2)%3];
00418 
00419             std::deque<vtkEDFEdge>::iterator beforeit;
00420             if(contourit == contour.begin())
00421               beforeit = contour.end()-1;
00422             else
00423               beforeit = contourit - 1;
00424             vtkEDFEdge& beforeedge = *beforeit;
00425 
00426             std::deque<vtkEDFEdge>::iterator afterit;
00427             if(contourit == contour.end()-1)
00428               afterit = contour.begin();
00429             else
00430               afterit = contourit + 1;
00431             vtkEDFEdge& afteredge = *afterit;
00432 
00433             if(beforeedge == tbeforeedge)
00434               {
00435               if(afteredge == tafteredge)
00436                 {// in this case, I am adding a triangle that is fully inside
00437                 // the contour. I need to remove the three edges from the
00438                 // contour.
00439                 if(contour.size() == 3)
00440                   {
00441                   contour.clear();
00442                   }
00443                 else
00444                   {
00445                   std::deque<vtkEDFEdge>::iterator lastit;
00446                   if(afterit == contour.end()-1)
00447                     lastit = contour.begin();
00448                   else
00449                     lastit = afterit + 1;
00450 
00451                   if(lastit < beforeit)
00452                     {
00453                     contour.erase(beforeit, contour.end());
00454                     contour.erase(contour.begin(), lastit);
00455                     }
00456                   else
00457                     {
00458                     contour.erase(beforeit, lastit);
00459                     }
00460                   }
00461                 }
00462               else
00463                 {// the edge before is the glued, remove the two adjacent edges
00464                 // and add the edge after
00465                 if(beforeit == contour.end()-1)
00466                   {
00467                   contour.erase(beforeit, contour.end());
00468                   contour.erase(contour.begin(), contour.begin()+1);
00469                   contour.push_back(tafteredge);
00470                   }
00471                 else
00472                   {
00473                   int index = beforeit - contour.begin();
00474                   contour.erase(beforeit, contourit+1);
00475                   contour.insert(contour.begin()+index, tafteredge);
00476                   }
00477                 }
00478               }
00479             else if(afteredge == tafteredge)
00480               {// the edge after is glued, removed the two glued edges and add
00481               // the edge new edge
00482               if(contourit == contour.end() -1)
00483                 {
00484                 contour.erase(contour.end() - 1);
00485                 contour.erase(contour.begin());
00486                 contour.push_back(tbeforeedge);
00487                 }
00488               else
00489                 {
00490                 int index = contourit - contour.begin();
00491                 contour.erase(contourit, afterit+1);
00492                 contour.insert(contour.begin()+index, tbeforeedge);
00493                 }
00494               }
00495             else
00496               {
00497               int index = contourit - contour.begin();
00498               contour.erase(contourit);
00499               contour.insert(contour.begin()+index, tbeforeedge);
00500               contour.insert(contour.begin()+index+1, tafteredge);
00501               }
00502             }
00503           vtkSmartPointer<vtkIdList> ids = vtkSmartPointer<vtkIdList>::New();
00504           std::deque<vtkEDFEdge>::iterator cit = contour.begin();
00505           while(cit != contour.end())
00506             {
00507             vtkEDFEdge& e = *cit;
00508             cit++;
00509             ids->InsertNextId(e.v0);
00510             }
00511 
00512           vtkIdType newCellId = outpd->InsertNextCell(VTK_POLYGON, ids);
00513           outpd->GetCellData()->CopyData(inpd->GetCellData(), cellId, newCellId);
00514           }
00515         }
00516       cellId = cellIdEnd - 1;
00517       }
00518     }
00519 }
00520 
00521 void  vtkEDFCutter::PrintSelf(ostream& os, vtkIndent indent)
00522 {
00523   this->Superclass::PrintSelf(os, indent);
00524 }
00525 
00526