Catmull–Clark subdivision surface: Difference between revisions
Catmull–Clark subdivision surface (view source)
Revision as of 12:15, 12 November 2023
, 6 months ago→{{header|Wren}}: Minor tidy
m (→{{header|Phix}}: made it p2js compatible) |
m (→{{header|Wren}}: Minor tidy) |
||
| (6 intermediate revisions by 4 users not shown) | |||
Line 1:
[[Category:Geometry]]
{{task|3D}}{{requires|Graphics}}
Implement the Catmull-Clark surface subdivision ([[wp:Catmull–Clark_subdivision_surface|description on Wikipedia]]), which is an algorithm that maps from a surface (described as a set of points and a set of polygons with vertices at those points) to another more refined surface. The resulting surface will always consist of a mesh of quadrilaterals.
Line 42 ⟶ 44:
For edges and vertices not next to a hole, the standard algorithm from above is used.
=={{header|C}}==
[[file:catmull-C.png|center]]
Only the subdivision part. The [[Catmull–Clark_subdivision_surface/C|full source]] is way too long to be shown here. Lots of macros, you'll have to see the full code to know what's what.
<
{
int i;
Line 138 ⟶ 139:
}
return nm;
}</
=={{header|Go}}==
{{trans|Python}}
Line 146:
See the original version for full comments.
<
import (
Line 431:
fmt.Printf("%2d\n", f)
}
}</
{{out}}
Line 487:
[ 4 23 13 20]
</pre>
=={{header|Haskell}}==
<
{-# LANGUAGE ScopedTypeVariables #-}
Line 799 ⟶ 798:
main :: IO ()
main = putStr . pShowSimpleMesh $ subdivide testCube</
{{out}}
<pre>Vertices:
Line 853 ⟶ 852:
(22,[24,15,5,8])
(23,[22,8,5,7])</pre>
=={{header|J}}==
<
havePoints=: e."1/~ i.@#
Line 878 ⟶ 876:
msh=. (,c0+mesh),.(,e0+edges i. meshEdges),.((#i.)~/$mesh),.,e0+_1|."1 edges i. meshEdges
msh;pts
)</
Example use:
<
points=: _1+2*#:i.8
mesh=: 1 A."1 I.(,1-|.)8&$@#&0 1">4 2 1
Line 914 ⟶ 912:
│ │ 0.555556 0.555556 _0.555556│
│ │ 0.555556 0.555556 0.555556│
└──────────┴─────────────────────────────┘</
=={{header|Julia}}==
<
# Point3f0 is a 3-tuple of 32-bit floats for 3-dimensional space, and all Points are 3D.
Line 1,062 ⟶ 1,059:
println("Press Enter to continue", readline())
</syntaxhighlight>
=={{header|Mathematica}}/{{header|Wolfram Language}}==
This implementation supports tris, quads, and higher polys, as well as surfaces with holes.
The function relies on three externally defined general functionality functions:
<
subSetQ[large_,small_List] := And@@(MemberQ[large,#]&/@small)
containing[groupList_,item_]:= Flatten[Position[groupList,group_/;subSetQ[group,item]]]
ReplaceFace[face_]:=Transpose[Prepend[Transpose[{#[[1]],face,#[[2]]}&/@Transpose[Partition[face,2,1,1]//{#,RotateRight[#]}&]],face]]</
subSetQ[small,large] is a boolean test for whether small is a subset of large. Note this is not a general purpose implimentation and only serves this purpose under the constrictions of the following program.
Line 1,083 ⟶ 1,079:
<
edges = DeleteDuplicates[Flatten[Partition[#,2,1,-1]&/@faces,1],Sort[#1]==Sort[#2]&];
avgFacePoints=Mean[Points[[#]]] &/@ faces;
Line 1,110 ⟶ 1,106:
newFaces = Flatten[Map[newIndex[#,{Points,edges,faces}]&,ReplaceFace/@faces,{-2}],1];
{newPoints,newFaces}
]</
The implimentation can be tested with polygons with and without holes by using the polydata
<
Function[iteration,
Graphics3D[(Polygon[iteration[[1]][[#]]]&/@iteration[[2]])]
]/@NestList[CatMullClark,{points,faces},3]//GraphicsRow</
[[File:CAM noholes 1.png]]
For a surface with holes the resulting iterative subdivision will be:
<
Function[iteration, Graphics3D[
(Polygon[iteration[[1]][[#]]] & /@ iteration[[2]])
]] /@ NestList[CatMullClark, {points, faces}, 3] // GraphicsRow</
[[File:CAM holes 1.png]]
This code was written in Mathematica 8.
=={{header|Nim}}==
{{trans|Go}}
{{trans|Python}}
<
import tables
Line 1,487 ⟶ 1,482:
s.add(fmt"{n:2d}")
s.add(']')
echo s</
{{out}}
Line 1,542 ⟶ 1,537:
[ 2, 20, 13, 17]
[ 4, 23, 13, 20]</pre>
=={{header|OCaml}}==
{{incorrect|OCaml|wrong output data}}
Line 1,552 ⟶ 1,546:
In the [[Catmull–Clark subdivision surface/OCaml|sub-page]] there is also a program in OCaml+OpenGL which displays a cube subdivided 2 times with this algorithm.
<
let add3 (x1, y1, z1) (x2, y2, z2) (x3, y3, z3) =
Line 1,687 ⟶ 1,681:
(Dynar.to_array da_points,
Dynar.to_array new_faces)
;;</
=== Another implementation ===
Another implementation which should work with holes, but has only been tested on a cube
{{works with|OCaml|4.02+}}
<
let zero = { x = 0.0; y = 0.0; z = 0.0 }
let add a b = { x = a.x+.b.x; y = a.y+.b.y; z = a.z+.b.z }
Line 1,764 ⟶ 1,758:
Face [c 0 0 0; c 0 1 0; c 1 1 0; c 1 0 0]; Face [c 0 0 1; c 0 1 1; c 1 1 1; c 1 0 1] ] in
show_faces cube;
show_faces (catmull_clark cube)</
with output:
<pre>surface {
Line 1,800 ⟶ 1,794:
Face: (0.7500, 0.0000, 0.7500) (0.5556, -0.5556, 0.5556) (0.0000, -0.7500, 0.7500) (0.0000, 0.0000, 1.0000)
}</pre>
=={{header|Phix}}==
{{trans|Go}}
<!--<
<span style="color: #000080;font-style:italic;">-- demo\rosetta\Catmull_Clark_subdivision_surface.exw</span>
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>
Line 2,020 ⟶ 2,013:
<span style="color: #7060A8;">pp</span><span style="color: #0000FF;">(</span><span style="color: #000000;">outputPoints</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">pp</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">sq_sub</span><span style="color: #0000FF;">(</span><span style="color: #000000;">outputFaces</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">),{</span><span style="color: #004600;">pp_IntFmt</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"%2d"</span><span style="color: #0000FF;">})</span>
<!--</
{{out}}
<pre>
Line 2,074 ⟶ 2,067:
{ 4,23,13,20}}
</pre>
=={{header|Python}}==
<
"""
Line 2,633 ⟶ 2,625:
graph_output(output_points, output_faces)
</syntaxhighlight>
=={{header|Rust}}==
{{trans|Python}}
<syntaxhighlight lang="rust">
pub struct Vector3 {pub x: f64, pub y: f64, pub z: f64, pub w: f64}
pub struct Triangle {pub r: [usize; 3], pub(crate) col: [f32; 4], pub(crate) p: [Vector3; 3], n: Vector3, pub t: [Vector2; 3]}
pub struct Mesh{pub v: Vec<Vector3>, pub f: Vec<Triangle>}
impl Triangle{
pub fn new() -> Triangle {return Triangle {r: [0, 0, 0], col: [0.0; 4], p: [Vector3::new(0.0, 0.0, 0.0), Vector3::new(0.0, 0.0, 0.0), Vector3::new(0.0, 0.0, 0.0)], n: Vector3::new(0.0, 0.0, 0.0), t: [Vector2::new(0.0, 0.0), Vector2::new(0.0, 0.0), Vector2::new(0.0, 0.0)]}}
pub fn copy(&self) -> Triangle {return Triangle {r: self.r.clone(), col: self.col, p: [self.p[0].copy(), self.p[1].copy(), self.p[2].copy()], n: self.n.copy(), t: [self.t[0].copy(), self.t[1].copy(), self.t[2].copy()]}}
}
impl Vector3 {
pub fn new(x: f64, y: f64, z: f64) -> Vector3 {return Vector3 {x, y, z, w: 1.0}}
pub fn normalize(&mut self) {
let l = (self.x * self.x + self.y * self.y + self.z * self.z).sqrt();
self.x /= l;
self.y /= l;
self.z /= l;
}
pub fn dot_product(v1: &Vector3, v2: &Vector3) -> f64 {return v1.x * v2.x + v1.y * v2.y + v1.z * v2.z}
pub fn cross_product(v1: &Vector3, v2: &Vector3) -> Vector3 {return Vector3::new(v1.y * v2.z - v1.z * v2.y, v1.z * v2.x - v1.x * v2.z, v1.x * v2.y - v1.y * v2.x)}
pub fn intersect_plane(plane_n: &Vector3, plane_p: &Vector3, line_start: &Vector3, line_end: &Vector3, mut t: f64) -> Vector3 {
let mut p_n = plane_n.copy();
p_n.normalize();
let plane_d = -Vector3::dot_product(&p_n, plane_p);
let ad = Vector3::dot_product(line_start, &p_n);
let bd = Vector3::dot_product(line_end, &p_n);
t = (-plane_d - ad) / (bd - ad);
let line = line_end.copy() - line_start;
let line_i = line * t;
return line_start.copy() + &line_i;
}
pub fn copy(&self) -> Vector3 {return Vector3 {x: self.x, y: self.y, z: self.z, w: self.w}}
}
impl Mesh {
pub fn get_face_points(&self) -> Vec<Vector3> {
let mut face_points: Vec<Vector3> = Vec::new();
for curr_face in &self.f {
let mut face_point = Vector3::new(0.0, 0.0, 0.0);
for curr_point_index in curr_face.r {
let curr_point = &self.v[curr_point_index];
face_point += &curr_point
}
face_point /= curr_face.r.len() as f64;
face_points.push(face_point.copy());
}
return face_points;
}
pub fn get_edges_faces(&self) -> Vec<[f64; 7]> {
let mut edges: Vec<[usize; 3]> = Vec::new();
for face_num in 0..self.f.len() {
let face: Triangle = self.f[face_num].copy();
let num_points = face.p.len();
for point_index in 0..num_points {
let mut point_num_1 = 0;
let mut point_num_2 = 0;
if point_index < num_points - 1 {
point_num_1 = face.r[point_index];
point_num_2 = face.r[point_index + 1];
} else {
point_num_1 = face.r[point_index];
point_num_2 = face.r[0];
}
if point_num_1 > point_num_2 {
let temp = point_num_1;
point_num_1 = point_num_2;
point_num_2 = temp;
}
edges.push([point_num_1, point_num_2, face_num]);
}
}
edges.sort();
let num_edges = edges.len();
let mut index = 0;
let mut merged_edges: Vec<[f64; 4]> = Vec::new();
while index < num_edges {
let e1 = edges[index];
if index < num_edges - 1 {
let e2 = edges[index + 1];
if e1[0] == e2[0] && e1[1] == e2[1] {
merged_edges.push([e1[0] as f64, e1[1] as f64, e1[2] as f64, e2[2] as f64]);
index += 2;
} else {
merged_edges.push([e1[0] as f64, e1[1] as f64, e1[2] as f64, -1.0]);
index += 1;
}
} else {
merged_edges.push([e1[0] as f64, e1[1] as f64, e1[2] as f64, -1.0]);
index += 1
}
}
let mut edges_centers = Vec::new();
for me in merged_edges {
let p1 = self.v[me[0] as usize].copy();
let p2 = self.v[me[1] as usize].copy();
let cp: Vector3 = Mesh::center_point(&p1, &p2);
edges_centers.push([me[0] as f64, me[1] as f64, me[2] as f64, me[3] as f64, cp.x, cp.y, cp.z]);
}
return edges_centers;
}
pub fn get_edge_points(&self, edges_faces: &Vec<[f64; 7]>, face_points: &Vec<Vector3>) -> Vec<Vector3> {
let mut edge_points = Vec::new();
for edge in edges_faces {
let cp = Vector3::new(edge[4], edge[5], edge[6]);
let fp1: Vector3 = face_points[edge[2] as usize].copy();
let mut fp2: Vector3 = fp1.copy();
if edge[3] != -1.0 { fp2 = face_points[edge[3] as usize].copy() };
let cfp = Mesh::center_point(&fp1, &fp2);
let edge_point = Mesh::center_point(&cp, &cfp);
edge_points.push(edge_point);
}
return edge_points
}
pub fn get_avg_face_points(&self, face_points: &Vec<Vector3>) -> Vec<Vector3> {
let num_points = self.v.len();
let mut temp_points = Vec::new();
let mut div: Vec<i32> = Vec::new();
for _ in 0..num_points {
temp_points.push(Vector3::new(0.0, 0.0, 0.0));
div.push(0)
};
for face_num in 0..self.f.len() {
let mut fp = face_points[face_num].copy();
for point_num in self.f[face_num].r {
let tp = temp_points[point_num].copy();
temp_points[point_num] = tp + &fp;
div[point_num] += 1;
}
}
let mut avg_face_points: Vec<Vector3> = Vec::new();
for i in 0..temp_points.len() {
let tp: Vector3 = temp_points[i].copy();
let t = tp / (div[i] as f64);
avg_face_points.push(t.copy());
}
return avg_face_points;
}
pub fn get_avg_mid_edges(&self, edges_faces: &Vec<[f64; 7]>) -> Vec<Vector3> {
let num_points = self.v.len();
let mut temp_points = Vec::new();
let mut div: Vec<i32> = Vec::new();
for point_num in 0..num_points{ temp_points.push(Vector3::new(0.0, 0.0, 0.0)); div.push(0)}
for edge in edges_faces {
let cp = Vector3::new(edge[4], edge[5], edge[6]);
for point_num in [edge[0] as usize, edge[1] as usize] {
let tp = temp_points[point_num].copy();
temp_points[point_num] = tp + &cp;
div[point_num] += 1
}
}
let mut avg_mid_edges: Vec<Vector3> = Vec::new();
for i in 0..temp_points.len(){
let ame: Vector3 = temp_points[i].copy() / (div[i] as f64);
avg_mid_edges.push(ame)}
return avg_mid_edges
}
pub fn get_points_faces(&self) -> Vec<i32> {
let num_points = self.v.len();
let mut points_faces: Vec<i32> = Vec::new();
for point_num in 0..num_points{points_faces.push(0)}
for face_num in 0..self.f.len() {
for point_num in self.f[face_num].r {
points_faces[point_num] += 1;
}
}
return points_faces
}
pub fn get_new_points(&self, points_faces: &Vec<i32>, avg_face_points: &Vec<Vector3>, avg_mid_edges: &Vec<Vector3>) -> Vec<Vector3> {
let mut new_points: Vec<Vector3> = Vec::new();
for point_num in 0..self.v.len() {
let n = points_faces[point_num] as f64;
let m1 = (n - 3.0) / n;
let m2 = 1.0 / n;
let m3 = 2.0 / n;
let old_coords = self.v[point_num].copy();
let p1 = old_coords * m1;
let afp = avg_face_points[point_num].copy();
let p2 = afp * m2;
let ame = avg_mid_edges[point_num].copy();
let p3 = ame * m3;
let p4 = p1 + &p2;
let new_coords = p4 + &p3;
new_points.push(new_coords);
}
return new_points;
}
pub fn switch_nums(point_nums: [f64; 2]) -> [f64; 2] {
return if point_nums[0] < point_nums[1] { point_nums } else {[point_nums[1], point_nums[0]]}
}
pub fn get_key(points: [f64; 2]) -> String {
return points[0].to_string() + ";" + &*points[1].to_string();
}
pub fn subdivide(&mut self) {
let face_points = self.get_face_points();
let edges_faces = self.get_edges_faces();
let edge_points = self.get_edge_points(&edges_faces, &face_points);
let avg_face_points = self.get_avg_face_points(&face_points);
let avg_mid_edges = self.get_avg_mid_edges(&edges_faces);
let points_faces = self.get_points_faces();
let mut new_points = self.get_new_points(&points_faces, &avg_face_points, &avg_mid_edges);
let mut face_point_nums = Vec::new();
let mut next_point_num = new_points.len();
for face_point in face_points {
new_points.push(face_point);
face_point_nums.push(next_point_num);
next_point_num += 1;
}
let mut edge_point_nums: HashMap<String, usize> = HashMap::new();
for edge_num in 0..edges_faces.len() {
let point_num_1 = edges_faces[edge_num][0];
let point_num_2 = edges_faces[edge_num][1];
let edge_point = edge_points[edge_num].copy();
new_points.push(edge_point);
edge_point_nums.insert(Mesh::get_key([point_num_1, point_num_2]), next_point_num);
next_point_num += 1;
}
let mut new_faces = Vec::new();
for old_face_num in 0..self.f.len() {
let old_face = self.f[old_face_num].copy();
let a = old_face.r[0] as f64;
let b = old_face.r[1] as f64;
let c = old_face.r[2] as f64;
let face_point_abc = face_point_nums[old_face_num];
let edge_point_ab = *edge_point_nums.get(&*Mesh::get_key(Mesh::switch_nums([a, b]))).unwrap();
let edge_point_bc = *edge_point_nums.get(&*Mesh::get_key(Mesh::switch_nums([b, c]))).unwrap();
let edge_point_ca = *edge_point_nums.get(&*Mesh::get_key(Mesh::switch_nums([c, a]))).unwrap();
new_faces.push([a, edge_point_ab as f64, face_point_abc as f64, edge_point_ca as f64]);
new_faces.push([b, edge_point_bc as f64, face_point_abc as f64, edge_point_ab as f64]);
new_faces.push([c, edge_point_ca as f64, face_point_abc as f64, edge_point_bc as f64]);
}
self.f = Vec::new();
for face_num in 0..new_faces.len() {
let curr_face = new_faces[face_num];
let mut t1: Triangle = Triangle::new();
let mut t2: Triangle = Triangle::new();
t1.p[0] = Vector3::new(new_points[curr_face[0] as usize].x, new_points[curr_face[0] as usize].y, new_points[curr_face[0] as usize].z);
t1.p[1] = Vector3::new(new_points[curr_face[1] as usize].x, new_points[curr_face[1] as usize].y, new_points[curr_face[1] as usize].z);
t1.p[2] = Vector3::new(new_points[curr_face[2] as usize].x, new_points[curr_face[2] as usize].y, new_points[curr_face[2] as usize].z);
t2.p[0] = Vector3::new(new_points[curr_face[2] as usize].x, new_points[curr_face[2] as usize].y, new_points[curr_face[2] as usize].z);
t2.p[1] = Vector3::new(new_points[curr_face[3] as usize].x, new_points[curr_face[3] as usize].y, new_points[curr_face[3] as usize].z);
t2.p[2] = Vector3::new(new_points[curr_face[0] as usize].x, new_points[curr_face[0] as usize].y, new_points[curr_face[0] as usize].z);
t1.r = [curr_face[0] as usize, curr_face[1] as usize, curr_face[2] as usize];
t2.r = [curr_face[2] as usize, curr_face[3] as usize, curr_face[0] as usize];
self.f.push(t1);
self.f.push(t2);
}
self.v = new_points;
}
}
</syntaxhighlight>
=={{header|Tcl}}==
This code handles both holes and arbitrary polygons in the input data.
<
# Use math functions and operators as commands (Lisp-like).
Line 2,758 ⟶ 3,036:
list [concat $newPoints $facepoints $edgepoints] $newFaces
}</
The [[/Tcl Test Code|test code]] for this solution is available as well. The example there produces the following partial toroid output image:
<center>[[File:Tcl-Catmull.png]]</center>
=={{header|Wren}}==
{{trans|Go}}
{{libheader|Wren-dynamic}}
{{libheader|Wren-sort}}
{{libheader|Wren-math}}
{{libheader|Wren-fmt}}
<syntaxhighlight lang="wren">import "./dynamic" for Tuple, Struct
import "./sort" for Sort
import "./math" for Int
import "./fmt" for Fmt
var Point = Tuple.create("Point", ["x", "y", "z"])
var fields = [
"pn1", // point number 1
"pn2", // point number 2
"fn1", // face number 1
"fn2", // face number 2
"cp" // center point
]
var Edge = Tuple.create("Edge", fields)
var PointEx = Struct.create("PointEx", ["p", "n"])
var sumPoint = Fn.new { |p1, p2| Point.new(p1.x + p2.x, p1.y + p2.y, p1.z + p2.z) }
var mulPoint = Fn.new { |p, m| Point.new(p.x * m, p.y * m, p.z * m) }
var divPoint = Fn.new { |p, d| mulPoint.call(p, 1/d) }
var centerPoint = Fn.new { |p1, p2| divPoint.call(sumPoint.call(p1, p2), 2) }
var getFacePoints = Fn.new { |inputPoints, inputFaces|
var facePoints = List.filled(inputFaces.count, null)
var i = 0
for (currFace in inputFaces) {
var facePoint = Point.new(0, 0, 0)
for (cpi in currFace) {
var currPoint = inputPoints[cpi]
facePoint = sumPoint.call(facePoint, currPoint)
}
facePoints[i] = divPoint.call(facePoint, currFace.count)
i = i + 1
}
return facePoints
}
var getEdgesFaces = Fn.new { |inputPoints, inputFaces|
var edges = []
var faceNum = 0
for (face in inputFaces) {
var numPoints = face.count
for (pointIndex in 0...numPoints) {
var pointNum1 = face[pointIndex]
var pointNum2 = (pointIndex < numPoints-1) ? face[pointIndex+1] : face[0]
if (pointNum1 > pointNum2) {
var t = pointNum1
pointNum1 = pointNum2
pointNum2 = t
}
edges.add([pointNum1, pointNum2, faceNum])
}
faceNum = faceNum + 1
}
var cmp = Fn.new { |e1, e2|
if (e1[0] == e2[0]) {
if (e1[1] == e2[1]) return (e1[2] - e2[2]).sign
return (e1[1] - e2[1]).sign
}
return (e1[0] - e2[0]).sign
}
var numEdges = edges.count
Sort.quick(edges, 0, numEdges-1, cmp)
var eIndex = 0
var mergedEdges = []
while (eIndex < numEdges) {
var e1 = edges[eIndex]
if (eIndex < numEdges-1) {
var e2 = edges[eIndex+1]
if (e1[0] == e2[0] && e1[1] == e2[1]) {
mergedEdges.add([e1[0], e1[1], e1[2], e2[2]])
eIndex = eIndex + 2
} else {
mergedEdges.add([e1[0], e1[1], e1[2], -1])
eIndex = eIndex + 1
}
} else {
mergedEdges.add([e1[0], e1[1], e1[2], -1])
eIndex = eIndex + 1
}
}
var edgesCenters = []
for (me in mergedEdges) {
var p1 = inputPoints[me[0]]
var p2 = inputPoints[me[1]]
var cp = centerPoint.call(p1, p2)
edgesCenters.add(Edge.new(me[0], me[1], me[2], me[3], cp))
}
return edgesCenters
}
var getEdgePoints = Fn.new { |inputPoints, edgesFaces, facePoints|
var edgePoints = List.filled(edgesFaces.count, null)
var i = 0
for (edge in edgesFaces) {
var cp = edge.cp
var fp1 = facePoints[edge.fn1]
var fp2 = (edge.fn2 == -1) ? fp1 : facePoints[edge.fn2]
var cfp = centerPoint.call(fp1, fp2)
edgePoints[i] = centerPoint.call(cp, cfp)
i = i + 1
}
return edgePoints
}
var getAvgFacePoints = Fn.new { |inputPoints, inputFaces, facePoints|
var numPoints = inputPoints.count
var tempPoints = List.filled(numPoints, null)
for (i in 0...numPoints) tempPoints[i] = PointEx.new(Point.new(0, 0, 0), 0)
for (faceNum in 0...inputFaces.count) {
var fp = facePoints[faceNum]
for (pointNum in inputFaces[faceNum]) {
var tp = tempPoints[pointNum].p
tempPoints[pointNum].p = sumPoint.call(tp, fp)
tempPoints[pointNum].n = tempPoints[pointNum].n + 1
}
}
var avgFacePoints = List.filled(numPoints, null)
var i = 0
for (tp in tempPoints) {
avgFacePoints[i] = divPoint.call(tp.p, tp.n)
i = i + 1
}
return avgFacePoints
}
var getAvgMidEdges = Fn.new { |inputPoints, edgesFaces|
var numPoints = inputPoints.count
var tempPoints = List.filled(numPoints, null)
for (i in 0...numPoints) tempPoints[i] = PointEx.new(Point.new(0, 0, 0), 0)
for (edge in edgesFaces) {
var cp = edge.cp
for (pointNum in [edge.pn1, edge.pn2]) {
var tp = tempPoints[pointNum].p
tempPoints[pointNum].p = sumPoint.call(tp, cp)
tempPoints[pointNum].n = tempPoints[pointNum].n + 1
}
}
var avgMidEdges = List.filled(tempPoints.count, null)
var i = 0
for (tp in tempPoints) {
avgMidEdges[i] = divPoint.call(tp.p, tp.n)
i = i + 1
}
return avgMidEdges
}
var getPointsFaces = Fn.new { |inputPoints, inputFaces|
var numPoints = inputPoints.count
var pointsFaces = List.filled(numPoints, 0)
for (faceNum in 0...inputFaces.count) {
for (pointNum in inputFaces[faceNum]) {
pointsFaces[pointNum] = pointsFaces[pointNum] + 1
}
}
return pointsFaces
}
var getNewPoints = Fn.new { |inputPoints, pointsFaces, avgFacePoints, avgMidEdges|
var newPoints = List.filled(inputPoints.count, null)
for (pointNum in 0...inputPoints.count) {
var n = pointsFaces[pointNum]
var m1 = (n-3) / n
var m2 = 1 / n
var m3 = 2 / n
var oldCoords = inputPoints[pointNum]
var p1 = mulPoint.call(oldCoords, m1)
var afp = avgFacePoints[pointNum]
var p2 = mulPoint.call(afp, m2)
var ame = avgMidEdges[pointNum]
var p3 = mulPoint.call(ame, m3)
var p4 = sumPoint.call(p1, p2)
newPoints[pointNum] = sumPoint.call(p4, p3)
}
return newPoints
}
var switchNums = Fn.new { |pointNums|
if (pointNums[0] < pointNums[1]) return pointNums
return [pointNums[1], pointNums[0]]
}
var cmcSubdiv = Fn.new { |inputPoints, inputFaces|
var facePoints = getFacePoints.call(inputPoints, inputFaces)
var edgesFaces = getEdgesFaces.call(inputPoints, inputFaces)
var edgePoints = getEdgePoints.call(inputPoints, edgesFaces, facePoints)
var avgFacePoints = getAvgFacePoints.call(inputPoints, inputFaces, facePoints)
var avgMidEdges = getAvgMidEdges.call(inputPoints, edgesFaces)
var pointsFaces = getPointsFaces.call(inputPoints, inputFaces)
var newPoints = getNewPoints.call(inputPoints, pointsFaces, avgFacePoints, avgMidEdges)
var facePointNums = []
var nextPointNum = newPoints.count
for (facePoint in facePoints) {
newPoints.add(facePoint)
facePointNums.add(nextPointNum)
nextPointNum = nextPointNum + 1
}
var edgePointNums = {}
for (edgeNum in 0...edgesFaces.count) {
var pointNum1 = edgesFaces[edgeNum].pn1
var pointNum2 = edgesFaces[edgeNum].pn2
var edgePoint = edgePoints[edgeNum]
newPoints.add(edgePoint)
edgePointNums[Int.cantorPair(pointNum1, pointNum2)] = nextPointNum
nextPointNum = nextPointNum + 1
}
var newFaces = []
var oldFaceNum = 0
for (oldFace in inputFaces) {
if (oldFace.count == 4) {
var a = oldFace[0]
var b = oldFace[1]
var c = oldFace[2]
var d = oldFace[3]
var facePointAbcd = facePointNums[oldFaceNum]
var p = switchNums.call([a, b])
var edgePointAb = edgePointNums[Int.cantorPair(p[0], p[1])]
p = switchNums.call([d, a])
var edgePointDa = edgePointNums[Int.cantorPair(p[0], p[1])]
p = switchNums.call([b, c])
var edgePointBc = edgePointNums[Int.cantorPair(p[0], p[1])]
p = switchNums.call([c, d])
var edgePointCd = edgePointNums[Int.cantorPair(p[0], p[1])]
newFaces.add([a, edgePointAb, facePointAbcd, edgePointDa])
newFaces.add([b, edgePointBc, facePointAbcd, edgePointAb])
newFaces.add([c, edgePointCd, facePointAbcd, edgePointBc])
newFaces.add([d, edgePointDa, facePointAbcd, edgePointCd])
}
oldFaceNum = oldFaceNum + 1
}
return [newPoints, newFaces]
}
var inputPoints = [
Point.new(-1, 1, 1),
Point.new(-1, -1, 1),
Point.new( 1, -1, 1),
Point.new( 1, 1, 1),
Point.new( 1, -1, -1),
Point.new( 1, 1, -1),
Point.new(-1, -1, -1),
Point.new(-1, 1, -1)
]
var inputFaces = [
[0, 1, 2, 3],
[3, 2, 4, 5],
[5, 4, 6, 7],
[7, 0, 3, 5],
[7, 6, 1, 0],
[6, 1, 2, 4]
]
var outputPoints = inputPoints.toList
var outputFaces = inputFaces.toList
var iterations = 1
for (i in 0...iterations) {
var res = cmcSubdiv.call(outputPoints, outputFaces)
outputPoints = res[0]
outputFaces = res[1]
}
for (p in outputPoints) {
Fmt.aprint([p.x, p.y, p.z], 7, 4, "[]")
}
System.print()
for (f in outputFaces) {
Fmt.aprint(f, 2, 0, "[]")
}</syntaxhighlight>
{{out}}
<pre>
[-0.5556 0.5556 0.5556]
[-0.5556 -0.5556 0.5556]
[ 0.5556 -0.5556 0.5556]
[ 0.5556 0.5556 0.5556]
[ 0.5556 -0.5556 -0.5556]
[ 0.5556 0.5556 -0.5556]
[-0.5556 -0.5556 -0.5556]
[-0.5556 0.5556 -0.5556]
[ 0.0000 0.0000 1.0000]
[ 1.0000 0.0000 0.0000]
[ 0.0000 0.0000 -1.0000]
[ 0.0000 1.0000 0.0000]
[-1.0000 0.0000 0.0000]
[ 0.0000 -1.0000 0.0000]
[-0.7500 0.0000 0.7500]
[ 0.0000 0.7500 0.7500]
[-0.7500 0.7500 0.0000]
[ 0.0000 -0.7500 0.7500]
[-0.7500 -0.7500 0.0000]
[ 0.7500 0.0000 0.7500]
[ 0.7500 -0.7500 0.0000]
[ 0.7500 0.7500 0.0000]
[ 0.7500 0.0000 -0.7500]
[ 0.0000 -0.7500 -0.7500]
[ 0.0000 0.7500 -0.7500]
[-0.7500 0.0000 -0.7500]
[ 0 14 8 15]
[ 1 17 8 14]
[ 2 19 8 17]
[ 3 15 8 19]
[ 3 19 9 21]
[ 2 20 9 19]
[ 4 22 9 20]
[ 5 21 9 22]
[ 5 22 10 24]
[ 4 23 10 22]
[ 6 25 10 23]
[ 7 24 10 25]
[ 7 16 11 24]
[ 0 15 11 16]
[ 3 21 11 15]
[ 5 24 11 21]
[ 7 25 12 16]
[ 6 18 12 25]
[ 1 14 12 18]
[ 0 16 12 14]
[ 6 18 13 23]
[ 1 17 13 18]
[ 2 20 13 17]
[ 4 23 13 20]
</pre>
| |||