MeshKit  1.0
HarmonicMapper.cpp
Go to the documentation of this file.
00001 #include "HarmonicMapper.hpp"
00002 #include <iostream>
00003 #include <math.h>
00004 #include <map>
00005 
00006 namespace MeshKit {
00007 
00008 HarmonicMapper::HarmonicMapper(MKCore* core, vector<Vertex> &v, vector<Face> &t, vector<Edge> &e, vector<set<int> > a)
00009 {
00010   mk_core = core;
00011   vtx.insert(vtx.begin(), v.begin(), v.end());
00012   tri.insert(tri.begin(), t.begin(), t.end());
00013   edges.insert(edges.begin(), e.begin(), e.end());
00014   adj.insert(adj.begin(), a.begin(), a.end());
00015   
00016 }
00017 
00018 void HarmonicMapper::execute()
00019 {
00020         _iterative_map(0.001);
00021 
00022 }
00023 
00024 void HarmonicMapper::getUV(vector<Vertex> &v)
00025 {
00026         int count = -1;
00027         for (vector<Vertex>::iterator it = vtx.begin(); it != vtx.end(); it++){
00028                 count++;
00029                 if (it->onBoundary)  continue;
00030                 v[count].uv[0] = it->uv[0];
00031                 v[count].uv[1] = it->uv[1];
00032                 //std::cout << "HarmonicMapper index = " << count << "\tuv = {" << it->uv[0] << "," << it->uv[1] << "}\n";
00033         }       
00034 
00035 }
00036 
00037 //matrix method
00038 void HarmonicMapper::_map()
00039 {
00040         int n_interior = 0, n_boundary = 0;
00041         for (vector<Vertex>::iterator it = vtx.begin(); it != vtx.end(); it++){
00042                 if (it->onBoundary){
00043                         it->id = n_boundary;                    
00044                         n_boundary++;
00045                 }
00046                 else{
00047                         it->id = n_interior;
00048                         n_interior++;
00049                 }
00050         }
00051         
00052         //using Armadillo to solve linear equations
00053         mat A = zeros<mat>(n_interior, n_interior);
00054         mat B = zeros<mat>(n_interior, n_boundary);
00055         //setting up the matrix A
00056 
00057         for (vector<Vertex>::iterator it = vtx.begin(); it != vtx.end(); it++){
00058                 if (it->onBoundary) continue;
00059                 int index_i = it->id;
00060                 double sw = 0.0;
00061                 for (set<int>::iterator it_adj = adj[it->index].begin(); it_adj != adj[it->index].end(); it_adj++){
00062                         int index_vtx = -1;
00063                         if (edges[*it_adj].connect[0]->index == it->index)
00064                                 index_vtx = edges[*it_adj].connect[1]->index;
00065                         else
00066                                 index_vtx = edges[*it_adj].connect[0]->index;
00067                         int index_j = vtx[index_vtx].id;
00068                         double w = edges[*it_adj].e;
00069                         if (vtx[index_vtx].onBoundary)
00070                                 B(index_i, index_j) = w;
00071                         else
00072                                 A(index_i, index_j) = -1.0*w;
00073                         
00074                         sw += w;
00075                 }
00076                 A(index_i, index_i) = sw;
00077         }
00078         
00079         mat b = zeros<mat>(n_interior, 2);
00080         
00081         for (vector<Vertex>::iterator it = vtx.begin(); it != vtx.end(); it++){
00082                 if (!it->onBoundary)  continue;
00083                 b(it->id,0) = it->uv[0];
00084                 b(it->id,1) = it->uv[1];
00085         }
00086         mat c = B * b;//n_interior * 2
00087     mat uv_coord = solve(A, c);
00088         
00089         int count = -1;
00090         for (vector<Vertex>::iterator it = vtx.begin(); it != vtx.end(); it++){
00091                 if (it->onBoundary)  continue;
00092                 count++;
00093                 it->uv[0] = uv_coord(count, 0);
00094                 it->uv[1] = uv_coord(count, 1);
00095         }       
00096         
00097 
00098 }
00099 
00100 //iterative method
00101 void HarmonicMapper::_iterative_map(double epsilon)
00102 {
00103         //move interior vertices to its center of neighbors
00104         //boundary nodes are set to (0,0)
00105         vector<int> interior;
00106         for (std::vector<Vertex>::iterator it = vtx.begin(); it != vtx.end(); it++)
00107                 if (!it->onBoundary){
00108                         interior.push_back(it->index);
00109                         it->uv[0] = 0.0;
00110                         it->uv[1] = 0.0;
00111                 }
00112                                 
00113 
00114         while(true){
00115                 double error = -1.0e+10;                
00116                 
00117                 
00118                 for (std::vector<int>::iterator it = interior.begin(); it != interior.end(); it++){
00119                         Vector2D uv;
00120                         uv[0] = 0.0;
00121                         uv[1] = 0.0;
00122                         double weight = 0.0;
00123                         for (set<int>::iterator it_e = adj[*it].begin(); it_e != adj[*it].end(); it_e++){
00124                                 int adj_v = -1;
00125                                 if (edges[*it_e].connect[0]->index == (*it))
00126                                         adj_v = edges[*it_e].connect[1]->index;
00127                                 else
00128                                         adj_v = edges[*it_e].connect[0]->index;
00129                                 uv[0] += edges[*it_e].e*vtx[adj_v].uv[0];
00130                                 uv[1] += edges[*it_e].e*vtx[adj_v].uv[1];
00131                                 weight += edges[*it_e].e;
00132                         }
00133                         Vector2D pre_uv;
00134                         pre_uv[0] = vtx[*it].uv[0];
00135                         pre_uv[1] = vtx[*it].uv[1];
00136                         vtx[*it].uv[0] = uv[0]/weight;
00137                         vtx[*it].uv[1] = uv[1]/weight;
00138                         double v_err = sqrt(pow(pre_uv[0]-vtx[*it].uv[0],2)+pow(pre_uv[1]-vtx[*it].uv[1],2));
00139                         error = (v_err > error)? v_err : error;
00140                 }
00141                 
00142                 if (error < epsilon) break;
00143         }
00144 }
00145 
00146 
00147 HarmonicMapper::~HarmonicMapper()
00148 {
00149   std::cout << "It is over now in smoothing" << endl;
00150 }
00151 
00152 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines