Subversion Repositories gelsvn

Rev

Details | Last modification | View Log | RSS feed

Rev Author Line No. Line
688 khor 1
//
2
//  polygonize.cpp
3
//  PointReconstruction
4
//
5
//  Created by J. Andreas Bærentzen on 16/03/13.
6
//  Copyright (c) 2013 J. Andreas Bærentzen. All rights reserved.
7
//
8
 
9
#include "../Geometry/Implicit.h"
10
#include "polygonize.h"
11
#include "smooth.h"
12
#include "cleanup.h"
13
#include "triangulate.h"
14
 
15
using namespace std;
16
using namespace CGLA;
17
using namespace Geometry;
18
using namespace HMesh;
19
 
20
namespace {
21
    Vec3d xpf[] = {Vec3d(0,-0.5,-0.5),Vec3d(0,0.5,-0.5),Vec3d(0,0.5,0.5),Vec3d(0,-0.5,0.5)};
22
    Vec3d xmf[] = {Vec3d(0, 0.5,-0.5),Vec3d(0,-0.5,-0.5),Vec3d(0,-0.5,0.5),Vec3d(0,0.5,0.5)};
23
 
24
    Vec3d ypf[] = {Vec3d( 0.5,0, -0.5),Vec3d(-0.5,0, -0.5),Vec3d(-0.5,0, 0.5),Vec3d(0.5,0, 0.5)};
25
    Vec3d ymf[] = {Vec3d(-0.5,0, -0.5),Vec3d(0.5,0, -0.5),Vec3d(0.5,0, 0.5),Vec3d(-0.5,0, 0.5)};
26
 
27
    Vec3d zpf[] = {Vec3d(-0.5,-0.5,0),Vec3d(0.5,-0.5,0),Vec3d(0.5,0.5,0),Vec3d(-0.5,0.5,0)};
28
    Vec3d zmf[] = {Vec3d( 0.5,-0.5,0),Vec3d(-0.5,-0.5,0),Vec3d(-0.5,0.5,0),Vec3d(0.5,0.5,0)};
29
}
30
 
31
namespace HMesh
32
{
33
 
34
    void polygonize(const XForm& xform, const RGrid<float>& grid, std::vector<CGLA::Vec3d>& quad_vertices, float tau)
35
    {
36
        quad_vertices.clear();
37
        for(int i=0;i<xform.get_dims()[0];++i)
38
            for(int j=0;j<xform.get_dims()[1];++j)
39
                for(int k=0;k<xform.get_dims()[2];++k)
40
                {
41
                    Vec3i vox(i,j,k);
42
                    if(grid[vox] <= tau)
43
                    {
44
                        if(grid.in_domain(Vec3i(i+1,j,k)) && grid[Vec3i(i+1,j,k)] > tau)
45
                            for(int n=0;n<4;++n)
46
                                quad_vertices.push_back(xform.inverse(Vec3d(i+0.5, j,k) + xpf[n]));
47
                        if(grid.in_domain(Vec3i(i-1,j,k)) && grid[Vec3i(i-1,j,k)] > tau)
48
                            for(int n=0;n<4;++n)
49
                                quad_vertices.push_back(xform.inverse(Vec3d(i-0.5, j,k) + xmf[n]));
50
                        if(grid.in_domain(Vec3i(i,j+1,k)) && grid[Vec3i(i,j+1,k)] > tau)
51
                            for(int n=0;n<4;++n)
52
                                quad_vertices.push_back(xform.inverse(Vec3d(i, j+0.5,k) + ypf[n]));
53
                        if(grid.in_domain(Vec3i(i,j-1,k)) && grid[Vec3i(i,j-1,k)] > tau)
54
                            for(int n=0;n<4;++n)
55
                                quad_vertices.push_back(xform.inverse(Vec3d(i, j-0.5,k) + ymf[n]));
56
                        if(grid.in_domain(Vec3i(i,j,k+1)) && grid[Vec3i(i,j,k+1)] > tau)
57
                            for(int n=0;n<4;++n)
58
                                quad_vertices.push_back(xform.inverse(Vec3d(i, j,k+0.5) + zpf[n]));
59
                        if(grid.in_domain(Vec3i(i,j,k-1)) && grid[Vec3i(i,j,k-1)] > tau)
60
                            for(int n=0;n<4;++n)
61
                                quad_vertices.push_back(xform.inverse(Vec3d(i, j,k-0.5) + zmf[n]));
62
 
63
                    }
64
                }
65
 
66
    }
67
 
68
 
69
    void volume_polygonize(const XForm& xform, const Geometry::RGrid<float>& grid,
70
                           HMesh::Manifold& mani, float tau)
71
    {
72
        mani.clear();
73
        vector<Vec3d> quad_vertices;
74
        polygonize(xform, grid, quad_vertices, tau);
75
        vector<int> indices;
76
        vector<int> faces(quad_vertices.size()/4,4);
77
        for(int i=0;i<quad_vertices.size();++i)
78
            indices.push_back(i);
79
        mani.build(quad_vertices.size(),
80
                   quad_vertices[0].get(),
81
                   faces.size(),
82
                   &faces[0],
83
                   &indices[0]);
84
 
85
        stitch_more(mani, 1e-5);
86
        shortest_edge_triangulate(mani);
87
 
88
        float avg_edge_len=0;
89
        for(HalfEdgeIDIterator h = mani.halfedges_begin(); h != mani.halfedges_end();++h)
90
            avg_edge_len += length(mani, *h);
91
        avg_edge_len /= mani.no_halfedges();
92
        VolumetricImplicit imp(xform, grid);
93
        for(int iter=0;iter<4;++iter)
94
        {
95
            TAL_smoothing(mani, .25, 1);
96
            for(VertexIDIterator vid = mani.vertices_begin(); vid != mani.vertices_end(); ++vid)
97
                imp.push_to_surface(mani.pos(*vid),0,avg_edge_len*0.5);
98
        }
99
        mani.cleanup();
100
        cout << "Produced" << mani.no_faces() << " faces " << endl;
101
    }
102
 
103
}