Instructor : Dianne Hansford, PhDHomework 3 : Bezier

Patches and Surface Curvatures

The purpose of this assignment is to gain some experience with Bezier patches

and to understand surface curvature as a classifier of shape. In this nb, you will

work with Gaussian and mean curvature.Some code is provided. You must write

code where you see “TO DO”.

Build patches:

TO DO: Create four bicubic Bezier patch control net.

The first three patches, called elliptic, hyperbolic, parabolic, should define

patches that demonstrate these classifications of Gaussian curvature.

You can keep this process simple and define the x- and y-coordinates over a grid

in the xy-plane and then define the z-coordinates as a function of x and y to

achieve a desired shape.

The fourth patch should be a freeform design and called freeform. It must be

significantly different than the other surfaces and it should not be symmetric.

elliptic = {

};

hyperbolic = {

};

parabolic = {

};

freeform = {

};

The code in the following section loads the patches you defined into a Bezier

surface evaluator.We will study B-splines later, and you will learn that Bezier

methods are a special case of them.(Do not change this code.)In the code in the

rest of the nb, each of these will be referred to as the “surface evaluator”.

In[601]:=

surfelliptic = BSplineFunction[elliptic];

surfhyperbolic = BSplineFunction[hyperbolic];

surfparabolic = BSplineFunction[parabolic];

surffreeform = BSplineFunction[freeform];

Gaussian Curvature

Functions for computing Gaussian curvature.

To Do: Implement Gaussian curvature in this section. See “gaussK” below.

(* This global variable turns

on print statements when set to 1 *)

printflag = 0;

(* zero tolerance to avoid dividing by zero *)

zeroTol = 0.000000005;

(*

Compute the first partials du, dv, second partials duu, dvv, mixed partial (duv),

and the normal (nor) at a (u,v) position on the surface

*)

(* NO CHANGE TO THIS FUNCTION *)

partialsf_, u_, v_ :=  du = Derivative[1, 0][f][u, v];

dv = Derivative[0, 1][f][u, v];

2

duv = Derivative[1, 1][f][u, v];

duu = Derivative[2, 0][f][u, v];

dvv = Derivative[0, 2][f][u, v];

cp = Cross[du, dv];

ncp = Norm[cp];

Ifncp < zeroTol, ncp = zeroTol; nor = cp /ncp;

Ifprintflag ⩵ 1,  Print"Partial at u=", u, " v=", v; Print"du =", du, " du.du = ", du.du; Print"dv =", dv, " dv.dv = ", dv.dv; Print"duv =", duv, " du.dv = ", du.dv; Print"duu =", duu; Print"dvv =", dvv; Print"nor =", nor

; Returndu, dv, duu, dvv, duv, nor

;

(*

Compute the first and second fundamental forms at a point

on the surface correspoinding to (u,v) in the domain.

Input the derivative information in this order: du,

dv, duu, dvv, duv, nor

where du is the first partial in the u-direction and

duu is the second partial in the v-direction.

3

nor is the normal vector at the point on the surface.

*)

(* NO CHANGE TO THIS FUNCTION *)

fundamentalForms[derivatives_] := 

(* assign the input to easier to use variables *)

du = derivatives〚1〛; dv = derivatives〚2〛; duu = derivatives〚3〛; dvv = derivatives〚4〛; duv = derivatives〚5〛; nor = derivatives〚6〛;

firstFundForm = Detdu.du, du.dv, du.dv, dv.dv; secondFundForm =

Detnor.duu, nor.duv, nor.duv, nor.dvv;

Ifprintflag ⩵ 1,  Print"First fundamental form =", firstFundForm; Print"Second fundamental form =", secondFundForm

; ReturnfirstFundForm, secondFundForm

;

(*

Compute the Gaussian curvature at a point on the surface.

Input the surface evaluator and the (u,v) pair

*)

(* TO DO: WRITE THIS FUNCTION *)

4

gaussK[f_, u_, v_] := 

(* Check for division by zero. If divisor is small,

set it to the zero tolerance.

Exchange "xxx" for the divisor *)

(*

IfAbs[xxx]< zeroTol,

Print"Warning! divisor is small"; xxx = zeroTol;

*)

Return(* gaussian curvature *) ;

(*

Calculate minimum and maximum computed Gaussian

curvatures over the surface for use in the color map.

Input the surface evaluator.

*)

(* NO CHANGE TO THIS FUNCTION *)

setminmax[f_] := 

(* Extract minimum and maximum computed Gaussian

curvatures over the surface for domain [0,1] x [0,1]

*) mingc = NMinValue

gaussK[f, uu, vv], 0 ≤ uu ≤ 1, 0 ≤ vv ≤ 1, {uu, vv}; maxgc = NMaxValue

5

gaussK[f, uu, vv], 0 ≤ uu ≤ 1, 0 ≤ vv ≤ 1, {uu, vv};

(* Print"min and max Gaussian curvature: n mingc = ",

mingc," maxgc = ",maxgc; *)

(*

Using a zero-tolerance,

check if mingc equals maxgc. If so,

modify the maxgc by the tolerance.

This will allow the color

map definition to work in all cases.

*)

zeroTol = 0.000005;

Ifmaxgc – mingc < zeroTol,

Print"Setting maxgc to +zeroTol"; maxgc = mingc + zeroTol;

Returnmingc, maxgc ;

Plot the surfaces with Gaussian curvature as a texture map.

To Do: Choose a different color scheme, and use it for each of the four plots.

It must be a diverging color scheme (like TemperatureMap)

To Do: Under each plot, create a formatted text cell and in this cell, describe the

patch (shape and control points) and the Gaussian curvature, making an

argument as to why this plot is correct. Consider commenting on the color map.

Plot the elliptic surface and its Gaussian curvature

6

In[611]:=

(* Calculate the min and max Gaussian curvature *)

{mingc, maxgc} = setminmax[surfelliptic];

(* Plot the surface with Gaussian curvature texture map *)

(* Rescale is used to map the Gaussian

curvature into the color ramp in [0,1] *)

plotelliptic = ParametricPlot3Dsurfelliptic[u, v], {u, 0, 1},

{v, 0, 1}, Mesh → False, PlotRange → All, ColorFunction → Function[{x, y, z, u, v}, ColorData["TemperatureMap"][

Rescale[gaussK[surfelliptic, u, v], {mingc, maxgc}]]],

PlotLegends → BarLegend[{"TemperatureMap", {mingc, maxgc}},

LegendLabel → "K"];

(* Show the surface with Gaussian curvature

and the control net with control points *)

gaussEllipticPlot = Show[plotelliptic,

Graphics3D[{PointSize[Large], Black, Map[Point, elliptic]}],

Graphics3D[{Thick, Gray, Line[elliptic],

Gray, Line[Transpose[elliptic]]}]]

Print Style"Color map of Gaussian curvature K", Large, Bold;

PrintStyle"min K = ", Large, Bold, Style[NumberForm[mingc, 2], Large];

PrintStyle"max K = ", Large, Bold, Style[NumberForm[maxgc, 2], Large];

7

Your comments here

Plot the hyperbolic surface and its Gaussian curvature In[617]:=

{mingc, maxgc} = setminmax[surfhyperbolic];

plothyperbolic =

ParametricPlot3Dsurfhyperbolic[u, v], {u, 0, 1},

{v, 0, 1}, Mesh → False, PlotRange → All, ColorFunction → Function{x, y, z, u, v}, ColorData["TemperatureMap"]

RescalegaussK[surfhyperbolic, u, v], mingc, maxgc, PlotLegends → BarLegend[{"TemperatureMap", {mingc, maxgc}},

LegendLabel → "K"];

gaussHyperbolicPlot = Show[plothyperbolic,

Graphics3D[{PointSize[Large], Black, Map[Point, hyperbolic]}],

Graphics3D[{Thick, Gray, Line[hyperbolic],

Gray, Line[Transpose[hyperbolic]]}]]

Print Style"Color map of Gaussian curvature K", Large, Bold;

PrintStyle"min K= ", Large, Bold, Style[NumberForm[mingc, 2], Large, Bold];

PrintStyle"max K= ", Large, Bold, Style[NumberForm[maxgc, 2], Large, Bold];

Your comments here

8

Plot the parabolic surface and its Gaussian curvature In[623]:=

{mingc, maxgc} = setminmax[surfparabolic];

plotparabolic =

ParametricPlot3Dsurfparabolic[u, v], {u, 0, 1}, {v, 0, 1},

Mesh → False, PlotRange → All, ColorFunction → Function{x, y, z, u, v}, ColorData["TemperatureMap"]

RescalegaussK[surfparabolic, u, v], mingc, maxgc, PlotLegends → BarLegend[{"TemperatureMap", {mingc, maxgc}},

LegendLabel → "K"];

gaussParabolicPlot = Show[plotparabolic,

Graphics3D[{PointSize[Large], Black, Map[Point, parabolic]}],

Graphics3D[{Thick, Gray, Line[parabolic],

Gray, Line[Transpose[parabolic]]}]]

Print Style"Color map of Gaussian curvature K", Large, Bold;

PrintStyle"min K = ", Large, Bold, Style[NumberForm[mingc, 2], Large, Bold];

PrintStyle"max K = ", Large, Bold, Style[NumberForm[maxgc, 2], Large, Bold];

Your comments here

Plot the freeform surface and its Gaussian curvature

9

In[629]:=

{mingc, maxgc} = setminmax[surffreeform];

plotfreeform =

ParametricPlot3Dsurffreeform[u, v], {u, 0, 1}, {v, 0, 1},

Mesh → False, PlotRange → All, ColorFunction → Function{x, y, z, u, v}, ColorData["TemperatureMap"]

RescalegaussK[surffreeform, u, v], mingc, maxgc, PlotLegends → BarLegend[{"TemperatureMap", {mingc, maxgc}},

LegendLabel → "K"];

gaussFreeFormPlot = Show[plotfreeform,

Graphics3D[{PointSize[Large], Black, Map[Point, freeform]}],

Graphics3D[{Thick, Gray, Line[freeform],

Gray, Line[Transpose[freeform]]}]]

Print Style"Color map of Gaussian curvature K", Large, Bold;

PrintStyle"min K = ", Large, Bold, Style[NumberForm[mingc, 2], Large, Bold];

PrintStyle"max K = ", Large, Bold, Style[NumberForm[maxgc, 2], Large, Bold];

Your comments here

Mean Curvature

To Do: Implement mean curvature

(*

Compute the mean curvature at a point on the surface.

Input the surface evaluator and the (u,v) pair

10

*)

(* WRITE THIS FUNCTION *)

meanCurvature[f_, u_, v_] := 

Return(* mean curvature *) ;

(*

Calculate minimum and maximum computed mean

curvatures over the surface for use in the color map.

Input the surface evaluator.

*)

(* NO CHANGE TO THIS FUNCTION *)

setminmaxMean[f_] := 

(* Extract minimum and maximum computed Gaussian

curvatures over the surface for domain [0,1] x [0,1]

*) minmc = NMinValuemeanCurvature[f, uu, vv],

0 ≤ uu ≤ 1, 0 ≤ vv ≤ 1, {uu, vv}; maxmc = NMaxValuemeanCurvature[f, uu, vv],

0 ≤ uu ≤ 1, 0 ≤ vv ≤ 1, {uu, vv};

(*Print"min and max mean curvature: n minmc = ",

minmc," maxmc = ",maxmc; *)

(*

11

Using a zero-tolerance,

check if minmc equals maxmc. If so,

modify the maxmc by the tolerance.

This will allow the color

map definition to work in all cases.

*)

eps = 0.000005;

Ifmaxmc – minmc < eps,

Print"Setting maxmc to +eps"; maxmc = minmc + eps;

Returnminmc, maxmc ;

To Do: Repeat the visualizations of the four surface with mean curvature as the

texture map.

Use the same color scheme that you chose above.

To Do: Under each plot, create a formatted text cell and in this cell, describe the

patch and the mean curvature, making an argument as to why this plot is correct.

Plot the elliptic surface and its mean curvature

Plot the hyperbolic surface and its mean curvature

Plot the parabolic surface and its mean curvature

Plot the freeform surface and its mean curvature

12

Use GraphicsRow to examine Gaussian and mean curvature plots side-by-

side.To Do: Below each plot in a formatted cell, briefly compare the plots

.Does one curvature measure reveal more about the shape than the other?

In[661]:=

GraphicsRowgaussEllipticPlot, meanEllipticPlot In[662]:=

GraphicsRowgaussHyperbolicPlot, meanHyperbolicPlot In[663]:=

GraphicsRowgaussParabolicPlot, meanParabolicPlot In[664]:=

GraphicsRowgaussFreeFormPlot, meanFreeFormPlot

13

In[632]:=

ClearAllEvaluate$Context <> "*"

Basic Definition In[633]:=

colorInterpolation[startColor_, endColor_, numCurves_, currentIndex_] :=

Moduler1, g1, b1, r2, g2, b2, fraction, {r1, g1, b1} = startColor;

{r2, g2, b2} = endColor;

(*Calculate the fraction representing the current position between the curves*)

fraction = (currentIndex – 1)  (numCurves – 1);

(*Interpolate the color components separately*)

RGBColorr1 + fraction (r2 – r1), g1 + fraction (g2 – g1), b1 + fraction (b2 – b1)

(*Example usage*)

startColor = {1, 0, 0}; (*Red*)

endColor = {0, 1, 0}; (*Blue*)

numCurves = 9; (*Number of curves in the design*)

(*Calculate colors for three curve pieces*)

curveColors = TablecolorInterpolation[startColor, endColor, numCurves, i], i, numCurves Out[637]=

{ , , , , , , , , }

In[638]:=

globalToLocal[a_, b_, u_] := (u – a)  (b – a)

In[639]:=

bezierPoints = 

{0, 0}, {1, 0} , {2, 0}, {3, 0},

{3, 0}, 4, -5, 5, 5, 6, 0,

{6, 0}, {7, -5}, {8, -3}, 9, 0,

{9, 0}, 10, 3, {11, 0}, 12, 0,

{12, 0}, {14, 0}, {11, 3}, {13, 2},

{{13, 2}, {12, 1}, {15, 1}, {13, 2}},

{14, 3}, 15, 4, 16, 4, 17, 0,

{17, 0}, 18, – 10, -5, 0, {0, 0}

;

Below is very important. It’s the core part: a function that can generate the bezire

curve.

In[640]:=

f[u_] := ModuleindexInbezierPoints = Floor[u] + 1, localVar = u – Floor[u], pts,

pts = bezierPointsindexInbezierPoints; Sumptsi + 1 × BernsteinBasis[3, i, localVar], i, 0, 3



In[641]:=

g1 = ParametricPlotf[t], t, 0, Length@bezierPoints

Out[641]=

5 10 15

-4

-2

2

In[642]:=

bezierLinesGraph =

TableGraphicscurveColorsi, LinebezierPointsi, i, Length[bezierPoints];

In[643]:=

controlPoints = TableGraphicsPointSize[0.02], curveColorsi, PointbezierPointsi,

i, Length@bezierPoints;

In[644]:=

g2 = ParametricPlotf[t], t, 0, Length@bezierPoints, ColorFunction → Function[{x, y, u}, curveColors〚Floor[u] + 1〛], ColorFunctionScaling → False;

2

Plot1

In[645]:=

Show[bezierLinesGraph, g2, controlPoints]

Out[645]=

Plot2

In[646]:=

g1

Out[646]=

5 10 15

-4

-2

2

Plot3 (wait)

In[647]:=

initRec = {-1, 1}, {1, -1};

3

In[648]:=

GraphicsBlue, Disk [{0, 0}, 1.5](*添加蓝⾊填充的Disk*),

Red, Rectangle@@ initRec, Black, Arrow[{{0, 0}, {1, 0}}], Green,

Arrow{0, 0}, {0, 1}, (*添加垂直箭头*)PointSize[0.02], Point[{0, 0}]

Out[648]=

This is very important: The derivative of bazier curve. Ref:

In[649]:=

"https://computergraphics.stackexchange.com/questions/10551/how-to-take-the-derivative-of-

a-b%C3%A9zier-curve#:~:text=We%20obtain%20a%20quadratic%20B,control%20points%20

scaled%20by%20n.";

In[650]:=

firstD[u_] :=

ModuleindexInbezierPoints = Floor[u] + 1, localVar = u – Floor[u], diffs, two,

diffs = Differences /@ bezierPointsindexInbezierPoints; two = 3 Sumdiffsi + 1 × BernsteinBasis2, i, localVar, i, 0, 2;

(*Don't need to convert to a number*)

two〚1〛, two〚2〛



In[651]:=

firstD[0.5]

Out[651]=

{3., 0.}

4

In[652]:=

SecondD[u_] :=

ModuleindexInbezierPoints = Floor[u] + 1, localVar = u – Floor[u], diffs, diffs2, two,

diffs = Differences /@ bezierPointsindexInbezierPoints; diffs2 = Differences /@ diffsindexInbezierPoints; two = 6 Sumdiffsi + 1 × BernsteinBasis1, i, localVar, i, 0, 1;

(*Don't need to convert to a number*)

two〚1〛, two〚2〛



Plot4

In[653]:=

ManipulateModulecurve, tangent, tangentVector, u, perpendicularVector,

translationMatrix, rotationMatrix, movingRectangle, u = uValue;

curve = f[u];

tangent = firstD[u];

tangentVector = tangent;

perpendicularVector = {-tangentVector〚2〛, tangentVector〚1〛}; translationMatrix = TranslationTransform[curve];

rotationMatrix = RotationTransform[{{1, 0}, tangentVector}];

movingRectangle =

GeometricTransformation[Rectangle@@ initRec, translationMatrix.rotationMatrix];

blueDisk = Disk[curve, 1.5];

ShowGraphicsBlue, blueDisk (*添加蓝⾊填充的Disk*), Red,

movingRectangle, Black, Arrowcurve, curve + Normalize[tangentVector], Green, Arrowcurve, curve + Normalize[perpendicularVector], (*添加垂直箭头*)PointSize[0.02], Point[curve], g1,

uValue, 0, "Parameter u", 0, Length@bezierPoints – 0.02,

0.001, Appearance → "Labeled"

Out[653]=

Parameter u 0

Show , g1

5

In[1116]:=

Clear["Global`*"]

CSE 477

Instructor: Dianne Hansford

Homework 2:

Exploring Least Squares Approximation with Bezier Curves

The purpose of this homework is to explore least squares approximation with Bezier curves. The

program flow is as follows.1) Input data points to be approximated 2) Uniform or chord length

parameters are created for the data.3) Input the degree n of the approximation4) Approximate the

data with a degree n Bezier curve

This program finds the approximation with both sets of parameter values so they may be

displayed together.5) Measure error between the Bezier curve and the input data.6) Display the

input, approximation, errors7) Document your observations. Details for each of these steps are

given below in the blue sections. You must complete the code when instructed to do so. I have

supplied some code and graphics for you.See the Observations section for instructions for the final

report: running the data sets and reporting.Instructions on how to get started will be given in class.

Guidelines:

Turn in via Canvas

Name your file lastname_HW2.nb

Below each problem description, add your MM code

Comment your notebook

Work independently

Bezier Curve Evaluator

(You do not need to change this routine.)

In[1117]:=

curve[bpts_, deg_, t_] := SumBernsteinBasis[deg, i, t] * bptsi + 1, i, 0, deg;

Data Sets

Below are some sample data sets

Create numdata = number of data to be approximated

Create data (x,y) or (x,y,z) points to be approximated

Load one data set into “data”

In[1118]:=

(* uniformly sampled line *)

bpts = {0, 0}, {1, 1}, {2, 2}, {3, 3};

dataUniformline = Tablecurvebpts, Length[bpts] – 1, t, t, 0, 1, 1.0 / 10.0;

(* Noisy line *)

dataNoisyLine = Tablet, t + 0.2 * RandomReal[], t, 0, 1, 1.0 / 30.0;

(* nice curve: classic car fender sillouette *)

bpts = {0, 0}, {0, 1}, {0, 2}, {20, 3};

dataGF = Tablecurvebpts, Length[bpts] – 1, t, t, 0, 1, 1.0 / 20.0;

(* This data set has a gap in the data —

suppose some part was occluded or data lost *)

dataGFgap = {0, 0}, 0.0025000000000000005`, 0.15000000000000002`, 0.020000000000000004`, 0.30000000000000004`, 0.06750000000000003`, 0.45000000000000007`, 0.16000000000000003`, 0.6000000000000002`,

0.3125`, 0.75`, 0.5400000000000003`, 0.9000000000000001`, 0.8575000000000002`, 1.05`, 1.2800000000000002`, 1.2000000000000002`, 1.8225000000000002`, 1.3500000000000003`, 5.4925000000000015`, 1.95`, 6.860000000000001`, 2.1`, 8.4375`, 2.25`, 10.240000000000002`, 2.4000000000000004`, 12.282500000000002`, 2.5500000000000003`,

14.580000000000002`, 2.7`, 17.1475`, 2.8499999999999996`, {20, 3};

(* four points *)

bpts = {0, 0}, 0.0, 1.0, 2.0, 1.0, 2.0, 0.0;

dataFour = Tablecurvebpts, Length[bpts] – 1, t, t, 0, 1, 1.0 / 3.0;

(* This is a data set to simulate an outlier *)

dataOutlier = Tablet, Sin[t], t, 0, Pi, 1.0 / 10.0;

dataOutlier〚17, 2〛 = dataOutlier〚17, 2〛 – 0.5;

(* semici