laplacesmoother.cpp

   1 #include "laplacesmoother.h"
   2 #include <vtkCellLocator.h>
   3
   4 LaplaceSmoother::LaplaceSmoother()
   5 {
   6 }
   7
   8
   9 LaplaceSmoother::~LaplaceSmoother()
  10 {
  11 }
  12
  13 void LaplaceSmoother::operate()
  14 {
  15   DebugLevel=0;
  16   if(DebugLevel>10) cout<<"LaplaceSmoother reporting in."<<endl;
  17
  18   QVector<vtkIdType> AllCells;
  19   getAllSurfaceCells(AllCells, m_grid);
  20   QVector<vtkIdType> SelectedCells;
  21   getSurfaceCells(m_bcs, SelectedCells, m_grid);
  22
  23   EG_VTKDCC(vtkIntArray, cell_code, m_grid, "cell_code");
  24   createCellToCell(AllCells, c2c, m_grid);
  25
  26   QSet <vtkIdType> SelectedNodes;
  27   QSet <vtkIdType> InternalNodes;
  28   QSet <vtkIdType> ExternalNodes;
  29
  30   foreach(vtkIdType id_cell, SelectedCells)
  31   {
  32     vtkIdType N_pts, *pts;
  33     m_grid->GetCellPoints(id_cell, N_pts, pts);
  34     for(int i=0;i<N_pts;i++)
  35     {
  36       QSet <int> bc;
  37       foreach(vtkIdType C, n2c[pts[i]])
  38       {
  39         bc.insert(cell_code->GetValue(C));
  40       }
  41       if(DebugLevel>10) cout<<"pts[i]="<<pts[i]<<" and bc="<<bc<<endl;
  42       SelectedNodes.insert(pts[i]);
  43       if(bc.size()>1) ExternalNodes.insert(pts[i]);
  44       else
  45       {
  46         vtkIdType point=pts[i];
  47         QSet< int > NeighbourCells=n2c[point];
  48         vtkIdType start=*(NeighbourCells.begin());
  49         vtkIdType current=start;
  50         do
  51         {
  52           vtkIdType next=nextcell(current,point,c2c,m_grid);
  53           current=next;
  54         } while (current!=start && current!=-1);
  55         if(current==-1) ExternalNodes.insert(point);
  56         if(current==start) InternalNodes.insert(point);
  57       }
  58     }
  59   }
  60
  61   createNodeToNode(cells, nodes, _nodes, n2n, m_grid);
  62
  63   EG_VTKSP(vtkUnstructuredGrid,grid_tmp);
  64   EG_VTKSP(vtkUnstructuredGrid,m_grid_orig);
  65   makeCopy(m_grid, m_grid_orig);
  66
  67   double closestPoint[3];
  68   vtkIdType cellId;
  69   int subId;
  70   double dist2;
  71   vtkCellLocator* terminator=vtkCellLocator::New();
  72   terminator->SetDataSet(m_grid_orig);
  73   terminator->BuildLocator();
  74
  75   for(int i_iter=0;i_iter<NumberOfIterations;i_iter++)
  76   {
  77     if(DebugLevel>10) cout<<"i_iter="<<i_iter<<endl;
  78     makeCopy(m_grid, grid_tmp);
  79
  80     foreach(vtkIdType id_G,InternalNodes)
  81     {
  82       vec3_t G(0,0,0);
  83       foreach(int id_M,n2n[id_G])
  84       {
  85         vec3_t M;
  86         m_grid->GetPoint(id_M, M.data());
  87         G+=M;
  88       }
  89       G=(1./n2n[id_G].size())*G;
  90       vec3_t P;
  91       terminator->FindClosestPoint(G.data(),P.data(),cellId,subId,dist2);
  92       grid_tmp->GetPoints()->SetPoint(id_G, P.data());
  93     }
  94
  95     if(DebugLevel>10) cout << "SelectedNodes.size()=" << SelectedNodes.size() << endl;
  96     if(DebugLevel>10) cout << "InternalNodes.size()=" << InternalNodes.size() << endl;
  97     if(DebugLevel>10) cout << "ExternalNodes.size()=" << ExternalNodes.size() << endl;
  98     if(DebugLevel>10) cout << "InternalNodes=" << InternalNodes << endl;
  99
 100     makeCopy(grid_tmp,m_grid);
 101   }
 102   if(DebugLevel>10) cout_grid(cout,m_grid);
 103
 104 }