laplacesmoother.cpp

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