Subversion Repositories gelsvn

Rev

Rev 543 | Rev 572 | Go to most recent revision | Show entire file | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed

Rev 543 Rev 566
Line 67... Line 67...
67 
        double area_mixed = 0;
 67 
        double area_mixed = 0;
68 
        //For each triangle T from the 1-ring neighborhood of x
68 
        //For each triangle T from the 1-ring neighborhood of x
69 
        for(HalfEdgeWalker w = m.halfedgewalker(v); !w.full_circle(); w = w.circulate_vertex_cw()){
 69 
        for(HalfEdgeWalker w = m.halfedgewalker(v); !w.full_circle(); w = w.circulate_vertex_cw()){
70 
            double f_area = area(m, w.face());
 70 
            double f_area = area(m, w.face());
71 
 
 71 
 
- 
 
 72 
            Vec3d v0(m.pos(v));
72 
            Vec3d v1(m.pos(w.vertex()));
 73 
            Vec3d v1(m.pos(w.vertex()));
73 
            Vec3d v2(m.pos(w.next().vertex()));
 74 
            Vec3d v2(m.pos(w.next().vertex()));
74 
            Vec3d v0(m.pos(w.next().next().vertex()));
 - 
 
75 
 
 75 
 
76 
            double a0 = acos(dot(v1-v0, v2-v0)/(length(v1-v0)*length(v2-v0)));
 76 
            double a0 = acos(dot(v1-v0, v2-v0)/(length(v1-v0)*length(v2-v0)));
77 
            double a1 = acos(dot(v2-v1, v0-v1)/(length(v2-v1)*length(v0-v1)));
 77 
            double a1 = acos(dot(v2-v1, v0-v1)/(length(v2-v1)*length(v0-v1)));
78 
            double a2 = acos(dot(v0-v2, v1-v2)/(length(v0-v2)*length(v1-v2)));
 78 
            double a2 = acos(dot(v0-v2, v1-v2)/(length(v0-v2)*length(v1-v2)));
79 
 
 79 
 
Line 112... Line 112...
112 
            return;
 112 
            return;
113 
 
 113 
 
114 
        Vec3d vertex(m.pos(v));
 114 
        Vec3d vertex(m.pos(v));
115 
        curv_normal = Vec3d(0);
 115 
        curv_normal = Vec3d(0);
116 
        w_sum = 0;
 116 
        w_sum = 0;
117 
        for(HalfEdgeWalker walker = m.halfedgewalker(v); !walker.full_circle(); walker = walker.circulate_vertex_cw()){
 117 
        for(HalfEdgeWalker walker = m.halfedgewalker(v); !walker.full_circle(); walker = walker.circulate_vertex_ccw()){
118 
            Vec3d nbr(m.pos(walker.vertex()));
 118 
            Vec3d nbr(m.pos(walker.vertex()));
119 
            Vec3d left(m.pos(walker.next().vertex()));
 119 
            Vec3d left(m.pos(walker.next().vertex()));
120 
            Vec3d right(m.pos(walker.opp().prev().opp().vertex()));
 120 
            Vec3d right(m.pos(walker.opp().next().vertex()));
121 
 
 121 
 
122 
            double d_left = dot(normalize(nbr-left),
 122 
            double d_left = dot(cond_normalize(nbr-left),cond_normalize(vertex-left));
123 
                normalize(vertex-left));
 - 
 
124 
            double d_right = dot(normalize(nbr-right),
 123 
            double d_right = dot(cond_normalize(nbr-right),cond_normalize(vertex-right));
125 
                normalize(vertex-right));
 - 
 
126 
            double a_left  = acos(min(1.0, max(-1.0, d_left)));
 124 
            double a_left  = acos(min(1.0, max(-1.0, d_left)));
127 
            double a_right = acos(min(1.0, max(-1.0, d_right)));
 125 
            double a_right = acos(min(1.0, max(-1.0, d_right)));
128 
 
 126 
 
129 
            double w = 1.0/tan(a_left) + 1.0/tan(a_right);
 127 
            double w = 1.0/(1e-300+tan(a_left));
130 
 
 - 
 
- 
 
 128 
            w += 1.0/(1e-300+tan(a_right));
- 
 
 129 
//            double w = sin(a_left + a_right) / (1e-300 + sin(a_left)*sin(a_right));
131 
            curv_normal += w * (vertex-nbr);
 130 
            curv_normal += w * (nbr-vertex);
132 
            w_sum += w;
 131 
            w_sum += w;
133 
        }
 132 
        }
134 
 
 133 
 
135 
    }
 134 
    }
136 
 
 135 
 
Line 299... Line 298...
299 
    void mean_curvatures(const Manifold& m, VertexAttributeVector<double>& curvature, int smooth_steps)
 298 
    void mean_curvatures(const Manifold& m, VertexAttributeVector<double>& curvature, int smooth_steps)
300 
    {
 299 
    {
301 
        for(VertexIDIterator v = m.vertices_begin(); v != m.vertices_end(); ++v)
 300 
        for(VertexIDIterator v = m.vertices_begin(); v != m.vertices_end(); ++v)
302 
			if(!boundary(m,*v))
 301 
			if(!boundary(m,*v))
303 
			{
 302 
			{
304 
				Vec3d N = mean_curvature_normal(m, *v);
 303 
				Vec3d N = -mean_curvature_normal(m, *v);
305 
				curvature[*v] = length(N) * sign(dot(N,Vec3d(normal(m, *v))));
 304 
				curvature[*v] = length(N) * sign(dot(N,Vec3d(normal(m, *v))));
306 
			}
305 
			}
307 
        smooth_something_on_mesh(m, curvature, smooth_steps);
306 
        smooth_something_on_mesh(m, curvature, smooth_steps);
308 
    }
 307 
    }
309 
 
 308 
 
Line 330... Line 329...
330 
 
 329 
 
331 
            if(fabs(L[max_idx][max_idx])<fabs(L[min_idx][min_idx])) swap(max_idx, min_idx);
 330 
            if(fabs(L[max_idx][max_idx])<fabs(L[min_idx][min_idx])) swap(max_idx, min_idx);
332 
 
 331 
 
333 
            Mat3x3d frame_t = transpose(frame);
 332 
            Mat3x3d frame_t = transpose(frame);
334 
 
 333 
 
335 
            max_curv_direction[*v] = frame_t * Vec3d(Q[max_idx][0], Q[max_idx][1], 0);
 334 
            max_curv_direction[*v] = cond_normalize(frame_t * Vec3d(Q[max_idx][0], Q[max_idx][1], 0));
336 
 
 335 
 
337 
            min_curv_direction[*v] = frame_t * Vec3d(Q[min_idx][0], Q[min_idx][1], 0);
 336 
            min_curv_direction[*v] = cond_normalize(frame_t * Vec3d(Q[min_idx][0], Q[min_idx][1], 0));
338 
 
 337 
 
339 
            curvature[*v] = L[0][0]*L[1][1];
 338 
            curvature[*v] = L[0][0]*L[1][1];
340 
        }
 339 
        }
341 
    }
 340 
    }
342 
 
 341