Catmull–Clark subdivision surface: Difference between revisions

=={{header|Nim}}==

{{trans|Go}}

{{trans|Python}}

<lang Nim>import algorithm

import tables

const None = -1 # Index number used to indicate no data.

type

Point = array[3, float]

Face = seq[int]

Edge = object

pn1: int # Point number 1.

pn2: int # Point number 2.

fn1: int # Face number 1.

fn2: int # Face number 2.

cp: Point # Center point.

PointEx = object

p: Point

n: int

PointNums = tuple[pn1, pn2: int]

####################################################################################################

# Point operations.

func `+`(p1, p2: Point): Point =

## Adds points p1 and p2.

for i in 0..Point.high:

result[i] = p1[i] + p2[i]

#---------------------------------------------------------------------------------------------------

func `*`(p: Point; m: float): Point =

## Multiply point p by m.

for i in 0..Point.high:

result[i] = p[i] * m

#---------------------------------------------------------------------------------------------------

func `/`(p: Point; d: float): Point =

## Divide point p by d.

for i in 0..Point.high:

result[i] = p[i] / d

#---------------------------------------------------------------------------------------------------

func centerPoint(p1, p2: Point): Point =

## Return a point in the center of the segment ended by points p1 and p2.

for i in 0..Point.high:

result[i] = (p1[i] + p2[i]) / 2

####################################################################################################

func getFacePoints(inputPoints: seq[Point]; inputFaces: seq[Face]): seq[Point] =

## For each face, a face point is created which is the average of all the points of the face.

result.setLen(inputFaces.len)

for i, currFace in inputFaces:

var facePoint: Point

for currPointIndex in currFace:

# Add currPoint to facePoint. Will divide later.

facePoint = facePoint + inputPoints[currPointIndex]

result[i] = facePoint / currFace.len.toFloat

#---------------------------------------------------------------------------------------------------

func getEdgesFaces(inputPoints: seq[Point]; inputFaces: seq[Face]): seq[Edge] =

## Get list of edges and the one or two adjacent faces in a list.

## Also get center point of edge.

# Get edges for each face.

var edges: seq[array[3, int]]

for faceNum, face in inputFaces:

# Loop over index into face.

for pointIndex in 0..face.high:

# If not last point then edge is current point and next point

# and for last point, edge is current point and first point.

var pointNum1 = face[pointIndex]

var pointNum2 = if pointIndex < face.high: face[pointIndex + 1] else: face[0]

# Order points in edge by lowest point number.

if pointNum1 > pointNum2:

swap pointNum1, pointNum2

edges &= [pointNum1, pointNum2, faceNum]

# Sort edges by pointNum1, pointNum2, faceNum.

edges.sort(proc(a1, a2: array[3, int]): int =

result = cmp(a1[0], a2[0])

if result == 0:

result = cmp(a1[1], a2[1])

if result == 0:

result = cmp(a1[2], a2[2]))

# Merge edges with 2 adjacent faces:

# [pointNum1, pointNum2, faceNum1, faceNum2] or

# [pointNum1, pointNum2, faceNum1, None]

var eIndex = 0

var mergedEdges: seq[array[4, int]]

while eIndex < edges.len:

let e1 = edges[eIndex]

# Check if not last edge.

if eIndex < edges.high:

let e2 = edges[eIndex + 1]

if e1[0] == e2[0] and e1[1] == e2[1]:

mergedEdges &= [e1[0], e1[1], e1[2], e2[2]]

inc eIndex, 2

else:

mergedEdges &= [e1[0], e1[1], e1[2], None]

inc eIndex

else:

mergedEdges &= [e1[0], e1[1], e1[2], None]

inc eIndex

# Add edge centers.

for me in mergedEdges:

let p1 = inputPoints[me[0]]

let p2 = inputPoints[me[1]]

let cp = centerPoint(p1, p2)

result.add(Edge(pn1: me[0], pn2: me[1], fn1: me[2], fn2: me[3], cp: cp))

#---------------------------------------------------------------------------------------------------

func getEdgePoints(inputPoints: seq[Point];

edgesFaces: seq[Edge];

facePoints: seq[Point]): seq[Point] =

## For each edge, an edge point is created which is the average between the center of the

## edge and the center of the segment made with the face points of the two adjacent faces.

result.setLen(edgesFaces.len)

for i, edge in edgesFaces:

# Get center of two facepoints.

let fp1 = facePoints[edge.fn1]

# If not two faces, just use one facepoint (should not happen for solid like a cube).

let fp2 = if edge.fn2 == None: fp1 else: facePoints[edge.fn2]

let cfp = centerPoint(fp1, fp2)

# Get average between center of edge and center of facePoints.

result[i] = centerPoint(edge.cp, cfp)

#---------------------------------------------------------------------------------------------------

func getAvgFacePoints(inputPoints: seq[Point];

inputFaces: seq[Face];

facePoints: seq[Point]): seq[Point] =

## For each point calculate the average of the face points of the faces the

## point belongs to (avgFacePoints), create a list of lists of two numbers

## [facePointSum, numPoints] by going through the points in all the faces then

## create the avgFacePoints list of point by dividing pointSum(x, y, z) by numPoints

var tempPoints = newSeq[PointEx](inputPoints.len)

# Loop through faces, updating tempPoints.

for faceNum, pointNums in inputFaces:

let fp = facePoints[faceNum]

for pointNum in pointNums:

let tp = tempPoints[pointNum].p

tempPoints[pointNum].p = tp + fp

inc tempPoints[pointNum].n

# Divide to build the result.

result.setLen(inputPoints.len)

for i, tp in tempPoints:

result[i] = tp.p / tp.n.toFloat

#---------------------------------------------------------------------------------------------------

func getAvgMidEdges(inputPoints: seq[Point]; edgesFaces: seq[Edge]): seq[Point] =

## Return the average of the centers of edges the point belongs to (avgMidEdges).

## Create list with entry for each point. Each entry has two elements. One is a point

## that is the sum of the centers of the edges and the other is the number of edges.

## After going through all edges divide by number of edges.

var tempPoints = newSeq[PointEx](inputPoints.len)

# Go through edgesFaces using center updating each point.

for edge in edgesFaces:

for pointNum in [edge.pn1, edge.pn2]:

let tp = tempPoints[pointNum].p

tempPoints[pointNum].p = tp + edge.cp

inc tempPoints[pointNum].n

# Divide out number of points to get average.

result.setLen(inputPoints.len)

for i, tp in tempPoints:

result[i] = tp.p / tp.n.toFloat

#---------------------------------------------------------------------------------------------------

func getPointsFaces(inputPoints: seq[Point]; inputFaces: seq[Face]): seq[int] =

# Return the number of faces for each point.

result.setLen(inputPoints.len)

# Loop through faces updating the result.

for pointNums in inputFaces:

for pointNum in pointNums:

inc result[pointNum]

#---------------------------------------------------------------------------------------------------

func getNewPoints(inputPoints: seq[Point]; pointsFaces: seq[int];

avgFacePoints, avgMidEdges: seq[Point]): seq[Point] =

## m1 = (n - 3.0) / n

## m2 = 1.0 / n

## m3 = 2.0 / n

## newCoords = (m1 * oldCoords) + (m2 * avgFacePoints) + (m3 * avgMidEdges)

result.setLen(inputPoints.len)

for pointNum, oldCoords in inputPoints:

let n = pointsFaces[pointNum].toFloat

let (m1, m2, m3) = ((n - 3) / n, 1 / n, 2 / n)

let p1 = oldCoords * m1

let afp = avgFacePoints[pointNum]

let p2 = afp * m2

let ame = avgMidEdges[pointNum]

let p3 = ame * m3

let p4 = p1 + p2

result[pointNum] = p4 + p3

#---------------------------------------------------------------------------------------------------

func switchNums(pointNums: PointNums): PointNums =

## Return tuple of point numbers sorted least to most.

if pointNums.pn1 < pointNums.pn2: pointNums

else: (pointNums.pn2, pointNums.pn1)

#---------------------------------------------------------------------------------------------------

func cmcSubdiv(inputPoints: seq[Point];

inputFaces: seq[Face]): tuple[p: seq[Point], f: seq[Face]] =

## For each face, a face point is created which is the average of all the points of the face.

## Each entry in the returned list is a point (x, y, z).

let facePoints = getFacePoints(inputPoints, inputFaces)

let edgesFaces = getEdgesFaces(inputPoints, inputFaces)

let edgePoints = getEdgePoints(inputPoints, edgesFaces, facePoints)

let avgFacePoints = getAvgFacePoints(inputPoints, inputFaces, facePoints)

let avgMidEdges = getAvgMidEdges(inputPoints, edgesFaces)

let pointsFaces = getPointsFaces(inputPoints, inputFaces)

var newPoints = getNewPoints(inputPoints, pointsFaces, avgFacePoints, avgMidEdges)

#[Then each face is replaced by new faces made with the new points,

for a triangle face (a,b,c):

(a, edgePoint ab, facePoint abc, edgePoint ca)

(b, edgePoint bc, facePoint abc, edgePoint ab)

(c, edgePoint ca, facePoint abc, edgePoint bc)

for a quad face (a,b,c,d):

(a, edgePoint ab, facePoint abcd, edgePoint da)

(b, edgePoint bc, facePoint abcd, edgePoint ab)

(c, edgePoint cd, facePoint abcd, edgePoint bc)

(d, edgePoint da, facePoint abcd, edgePoint cd)

facePoints is a list indexed by face number so that is easy to get.

edgePoints is a list indexed by the edge number which is an index into edgesFaces.

Need to add facePoints and edgePoints to newPoints and get index into each.

Then create two new structures:

facePointNums: list indexes by faceNum whose value is the index into newPoints.

edgePointNums: dictionary with key (pointNum1, pointNum2) and value is index into newPoints.

]#

# Add face points to newPoints.

var facePointNums: seq[int]

var nextPointNum = newPoints.len # PointNum after next append to newPoints.

for facePoint in facePoints:

newPoints.add(facePoint)

facePointNums.add(nextPointNum)

inc nextPointNum

# Add edge points to newPoints.

var edgePointNums: Table[tuple[pn1, pn2: int], int]

for edgeNum, edgesFace in edgesFaces:

let pointNum1 = edgesFace.pn1

let pointNum2 = edgesFace.pn2

newPoints.add(edgePoints[edgeNum])

edgePointNums[(pointNum1, pointNum2)] = nextPointNum

inc nextPointNum

# newPoints now has the points to output. Need new faces.

#[Just doing this case for now:

for a quad face (a,b,c,d):

(a, edgePoint ab, facePoint abcd, edgePoint da)

(b, edgePoint bc, facePoint abcd, edgePoint ab)

(c, edgePoint cd, facePoint abcd, edgePoint bc)

(d, edgePoint da, facePoint abcd, edgePoint cd)

newFaces will be a list of lists where the elements are like this:

[pointNum1, pointNum2, pointNum3, pointNum4]

]#

var newFaces: seq[Face]

for oldFaceNum, oldFace in inputFaces:

# 4 points face.

if oldFace.len == 4:

let (a, b, c, d) = (oldFace[0], oldface[1], oldface[2], oldface[3])

let facePointAbcd = facePointNums[oldFaceNum]

let edgePointAb = edgePointNums[switchNums((a, b))]

let edgePointDa = edgePointNums[switchNums((d, a))]

let edgePointBc = edgePointNums[switchNums((b, c))]

let edgePointCd = edgePointNums[switchNums((c, d))]

newFaces &= @[a, edgePointAb, facePointAbcd, edgePointDa]

newFaces &= @[b, edgePointBc, facePointAbcd, edgePointAb]

newFaces &= @[c, edgePointCd, facePointAbcd, edgePointBc]

newFaces &= @[d, edgePointDa, facePointAbcd, edgePointCd]

result = (newPoints, newFaces)

#———————————————————————————————————————————————————————————————————————————————————————————————————

when isMainModule:

import strformat, strutils

let inputPoints = @[[-1.0, 1.0, 1.0],

[-1.0, -1.0, 1.0],

[1.0, -1.0, 1.0],

[1.0, 1.0, 1.0],

[1.0, -1.0, -1.0],

[1.0, 1.0, -1.0],

[-1.0, -1.0, -1.0],

[-1.0, 1.0, -1.0]]

let 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

var outputFaces = inputFaces

const Iterations = 1

for i in 1..Iterations:

(outputPoints, outputFaces) = cmcSubdiv(outputPoints, outputFaces)

for p in outputPoints:

echo fmt"[{p[0]: .4f}, {p[1]: .4f}, {p[2]: .4f}]"

echo ""

for nums in outputFaces:

var s = "["

for n in nums:

s.addSep(", ", 1)

s.add(fmt"{n:2d}")

s.add(']')

echo s</lang>

{{out}}

This is the same output as that of Go.

<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>