MagistrskodeloMentorica:izr.prof.dr.MarjetkaKnezLjubljana,2021 RAƒUNANJEUKRIVLJENOSTINAPROSTORSKIHTRIANGULACIJAH UNIVERZAVLJUBLJANIFAKULTETAZAMATEMATIKOINFIZIKOMatematika2.stopnjaPetraIvanaReberc

UNIVERZA V LJUBLJANI

FAKULTETA ZA MATEMATIKO IN FIZIKO Matematika 2. stopnja

Petra Ivana Reberc

RAƒUNANJE UKRIVLJENOSTI NA PROSTORSKIH TRIANGULACIJAH

Magistrsko delo

Mentorica: izr. prof. dr. Marjetka Knez

Ljubljana, 2021

Zahvala

Zahvaljujem se mami, Jeri, Katarini, Ja²otu, Maru²i, Niki, Adi in vsem, ki ste mi stali ob strani med leti ²tudija. Rada bi se zahvalila ²e mentorici prof. dr. Marjetki Knez za vso podporo, spodbudo in pomo£ pri nastajanju magistrske naloge.

Kazalo

Program dela vii

1 Uvod 1

2 Uvod v ukrivljenosti parametriziranih krivulj in ploskev 1 2.1 Fleksijska in torzijska ukrivljenost prostorskih parametriziranih krivulj 1

2.2 Parametrizirane ploskve in ukrivljenosti ploskev . . . 5

3 Ra£unanje ukrivljenosti s pomo£jo aproksimacije z Bézierjevimi ploskvami 15 3.1 Uvod v Bézierjeve krivulje in Bézierjeve ploskve . . . 15

3.2 Aproksimacija prostorske triangulacije z Bézierjevo ploskvijo po me- todi najmanj²ih kvadratov . . . 20

4 Ra£unanje ukrivljenosti z uporabo ploskev iz B-zlepkov 24 4.1 Aproksimacija prostorske triangulacije z zlepki . . . 30

5 Ra£unanje ukrivljenosti z uporabo metode prostorskih povpre£ij 32 5.1 Prostorsko povpre£je in 1-okolica ogli²£a . . . 33

5.2 Voronojeve in baricentri£ne celice . . . 34

5.3 Izra£un normale in povpre£ne ukrivljenosti . . . 38

5.4 Izra£un Gaussove ukrivljenosti . . . 43

5.5 Diskretni glavni ukrivljenosti . . . 51

5.6 Dolo£anje glavnih smeri . . . 52

Stvarno kazalo 55

Program dela

3D skeniranje je ena izmed najsodobnej²ih tehnologij, ki nam omogo£a prenos zi£ne oblike objekta v digitalen zapis, kar imenujemo digitalizacija. Z hitrim razvojem in dostopnostjo 3D skenerjev se v zadnjih letih pojavlja vse ve£ja potreba po nadaljnji analizi skeniranih objektov in avtomatskemu dolo£anju lastnosti le teh. Obi£ajno so izhodni podatki, ki jih dobimo iz digitalnih naprav, podani v obliki prostorskih triangulacij, pri dolo£anju lastnosti 3D ploskev, pa je ukrivljenost zagotovo ena od klju£nih lastnosti, ki jih je vredno analizirati. V magistrskem delu zato predstavite razli£ne pristope, tako globalne kot tudi lokalne narave, k aproksimaciji ukrivljenosti objektov, podanih v obliki prostorskih triangulacij. Pri tem si pomagajte s spodaj navedeno literaturo.

Osnovna literatura

[10] M. Meyer, M. Desbrun, P. Schröder in A. Barr, Discrete dierential-geometry operators for triangulated 2-manifolds, v: Visualization and mathematics III, Springer, 2003, str. 3557

[12] A. N. Pressley, Elementary dierential geometry, Springer Science & Business Media, 2010

Podpis mentorja:

Ra£unanje ukrivljenosti na prostorskih triangulacijah Povzetek

V delu so obravnavane razli£ne metode za ocenjevanje ukrivljenosti na dani pro- storski triangulaciji. Opisan je postopek aproksimacije triangulacije po metodi naj- manj²ih kvadratov z Bézierjevimi ploskvami iz tenzorskega produkta in z B-zlepki, na katerih lahko uporabimo standardne analiti£ne formule za izra£un ukrivljenosti.

Predstavljena je ²e metoda prostorskih povpre£ij, s katero izpeljemo ra£unsko nezah- tevne in elegantne formule za direkten izra£un ukrivljenosti v ogli²£ih triangulacije, pri £emer se upo²teva bliºnjo okolico izbranega ogli²£a.

Curvature computation on triangle meshes Abstract

The master thesis proposes various methods for estimating the curvature on tri- angle meshes in a 3D space. We describe the approximation of the triangle mesh with tensor product Bézier surfaces and B-spline surfaces, using the least squares method. We can then use standard formulas to compute curvature on the appro- ximating surface. Furthermore the method of spatial averages is presented, which oers computationally eective and elegant formulas for calculating the curvature at a given vertex in regard to its surroundings.

Math. Subj. Class. (2010): 65D07, 68U05

Klju£ne besede: ocena ukrivljenosti, prostorska triangulacija, metoda prostorskih povpre£ij, aproksimacija, Bézierjeva ploskev, ploskev iz B-zlepkov

Keywords: curvature estimation, triangle mesh, method of spatial averages, appro- ximation, Bézier surface, B-spline surface

1 Uvod

V delu se ukvarjamo z razli£nimi pristopi ra£unanja ukrivljenosti na prostorskih triangulacijah. Ko dobimo preko digitalnih naprav slike ploskev, so te pogosto po- dane v obliki prostorskih triangulacij. Te triangulacije so na£eloma samo zvezne, kar ni dovolj za uporabo formul za izra£un njihove ukrivljenosti, zato bi radi opisali moºne pristope, kako te ocene pridobimo. Ukrivljenost ploskve nas lahko zanima pri re²evanju razli£nih problemov. S spremljanjem sprememb ukrivljenosti ploskev na medicinskih slikah lahko bolj natan£no spremljamo potek bolezni, med seboj po- dobne izseke ploskve lahko identiciramo glede na ukrivljenost, itd. Kot bo kasneje opisano, se metode za ra£unanje aproksimacije ukrivljenosti v grobem delijo na dva dela. V prvem delu magistrske naloge se problem ra£unanja ukrivljenosti prevede na iskanje ustrezne aproksimacije ploskve po metodi najmanj²ih kvadratov. V tem sklopu bodo bolj podrobno predstavljene Bézierjeve ploskve iz tenzorskega produkta in ploskve izB-zlepkov. V drugem delu naloge pa je opisana t.i. metoda prostorskih povpre£ij, ki ponuja direktno oceno ukrivljenosti glede na okolico to£ke triangulacije, v kateri nas ukrivljenost zanima. Pri tem je bil osnovni vir £lanek [10]. Na samem za£etku pa se dotaknemo osnov krivulj, ploskev in ukrivljenosti, kar nam omogo£a bolj²e razumevanje problema. Prav tako je obravnavan eden klju£nih rezultatov di- ferencialne geometrije, Gauss-Bonnetov izrek, saj lahko z njegovo pomo£jo izpeljemo enostavne in ra£unsko nezahtevne formule za ocene ukrivljenosti.

2 Uvod v ukrivljenosti parametriziranih krivulj in ploskev

Razli£ne mere ukrivljenosti na ploskvah nam podajo bolj²e razumevanje lastnosti oblike ploskve. V tem razdelku so opisani osnovni pojmi iz diferencialne geometrije, na podlagi katerih so razvite metode za ocenjevanje ukrivljenosti na prostorski tri- angulaciji. V prvem razdelku deniramo torzijsko in eksijsko ukrivljenost krivulj, v drugem pa so opisane mere ukrivljenosti na ploskvah.

2.1 Fleksijska in torzijska ukrivljenost prostorskih parame- triziranih krivulj

Denicija 2.1. Denirajmo najprej prostorsko parametrizirano krivuljo. Naj bo I ⊂R odprt ali zaprt interval inr:I →R³ preslikava (vektorska funkcija)

r(t) = (x(t), y(t), z(t)).

Krivulja k je mnoºica to£k v zalogi vrednosti preslikaver.

Naj bo preslikava rpoljubnokrat odvedljiva na I, £e je ta odprt, sicer pa poljub- nokrat odvedljiva na vsakem odprtem intervalu I^′ ⊂I. ƒe velja ²e

ṙ (t)̸= 0,

kjer ṙ ozna£uje odvod po parametru t, je r regularna parametrizacija krivulje k. Parametrizirana krivulja je par, sestavljen iz mnoºice to£kkin preslikave r:I ⊂R³.

Opomba 2.2. V nadaljevanju bomo predpostavili, da imamo opravka z gladkimi krivuljami, izjeme bodo posebej poudarjene.

Opomba 2.3. Zaradi kraj²ega poimenovanja bomo dostikrat v nadaljevanju para- metrizaciji krivulje rekli kar krivulja.

Naj borparametrizacija krivulje. Zanima nas dolºina krivulje od to£ker(t₁) do to£ke r(t). Dolºina s je funkcija parametra t,

s(t) =

∫︂ t t1

√︁ṙ (τ)·ṙ (τ)dτ.

pri £emer · ozna£uje standardni skalarni produkt v R³. ƒe je parametrizacija r regularna, je s(t) strogo nara²£ajo£a funkcija, torej ima inverz na celotni domeni.

Tedaj lahko zapi²emo t kot funkcijo parametra s in reparametriziramo krivuljo z r(t(s)). Parameter s imenujemo naravni parameter. Predpostavimo naprej, da je krivuljarparametrizirana z naravnim parametrom. Posebnost te parametrizacije je, da veljaṙ·ṙ = 1. Vektorska funkcijaṙ v tem primeru dolo£a enotske vektorje v smeri tangente in jo ozna£imo s t. Vektorja t(s) inṫ(s) sta pravokotna za vsaks∈I, saj iz odvajanja enakosti t·t= 1sledi 2t·ṫ = 0. Smer, ki jo dolo£a ṫ,imenujemo smer glavne normale. Dolºino vektorja ṫ(s) imenujemo eksijska ukrivljenost krivulje v to£ki r(s). Ozna£imo jo z

κ= 1 ρ.

Videli bomo, da je eksijska ukrivljenost kroºnice z radijemρenaka 1/ρ. Velja torej κ = 1

ρ =√︁

ṫ·ṫ =√

¨r·r¨,

kjer so vsi odvodi po naravnem parametru. Smer, ki jo dolo£a vektorska funkcija b=t× ṫ

∥ṫ∥

imenujemo smer binormale. Opazimo, da je dolºina binormale v vsaki to£ki enaka1, saj sta si vektorja t(s)in _∥t^{ṫ (s)}_{̇ (s)∥} pravokotna in oba sta dolºine1. Ravnino, ki vsebuje dano to£ko, tangento in glavno normalo, imenujemo pritisnjena ravnina. Binor- mala je normala pritisnjene ravnine. Naj bo r(s) to£ka na krivulji pri parametrus. Tangentno premico na krivuljo v to£ki r(s)lahko opi²emo s parametrizacijo

R(u) =r(s) +ut(s), u∈R.

Pritisnjena ravnina, ki vsebuje r(s), je dolo£ena z implicitno ena£bo ((x, y, z)−r(s))·b= 0.

Denirajmo ²e normalno ravnino. To je ravnina, ki gre skozi dotikali²£e tangente in je pravokotna na tangento. Dolo£ena je z ena£bo

((x, y, z)−r(s))·t= 0.

Vemo ºe, da velja ∥b∥ = 1. Zato sta si vektorja b(s) in njegov odvod ḃ (s) po naravnem parametru pravokotna. Vpeljimo oznako

n= ṫ

∥ṫ∥.

Za vsak parameter s leºi vektor ḃ (s) v ravnini, ki je dolo£ena z vektorjema t(s) in n(s). Iz enakosti b=t×n sledi

ḃ =ṫ×n+t×ṅ

= 1

ρn×n+t×ṅ

=t×ṅ.

Vektorṅ (s)leºi v isti ravnini kotb(s)int(s), ker je pravokoten nan(s). Zato lahko zapi²emo n kot linearno kombinacijo binormale b in tangente t:

ṅ =αb+βt.

Sledi

t×ṅ =αt×b+βt×t=αt×b.

Iz lastnosti vektorskega produkta vemo, da ima t×b smer −n. Izpeljali smo ḃ =−ωn=−1

τn.

Skalarno funkcijo ω imenujemo torzijska ukrivljenost ali zvitost krivulje.

Za geometrijsko intuicijo si oglejmo, kaj se zgodi, ko sta torzijska in eksijska ukrivljenost krivulje enaki 0. Ko je eksijska ukrivljenost enaka 0, se tangenta ne spreminja, zato je krivulja premica. Ko je torzijska ukrivljenost enaka 0, pa je binormala konstantna, zato tangenta in glavna normala krivulje vedno leºita v isti ravnini, kar je moºno le, ko je krivulja ravninska. Z enostavnim razmislekom preverimo, da velja tudi obratno in ugotovitvi strnemo v slede£o trditev.

Trditev 2.4. Krivulja je premica natanko takrat, ko je njena ukrivljenost konstantno enaka 0. Krivulja je ravninska natanko takrat, ko je njena torzijska ukrivljenost enaka 0.

Denicija 2.5. Paroma ortonormirane vektorske funkcije(t,n,b)imenujemo spre- mljajo£i trieder krivulje oziroma Frenetovo ogrodje.

Ena£bi za eksijsko in torzijsko ukrivljenost smo zaenkrat izpeljali za primer, ko je krivulja parametrizirana z naravnim parametrom. Naj bo zdaj parametrizacija poljubna. Potem eksijsko ukrivljenost izra£unamo kot

κ= 1

ρ = ∥ṙ×r¨∥

∥ṙ∥³ , torzijsko ukrivljenost pa s formulo

1

τ = (ṙ ×r¨)· ...r

∥ṙ ר∥r ² .

Izpeljavo izpustimo, bralec pa jo lahko najde v [5, poglavje 3]. V tem viru je tudi pokazano, da je eksijska ukrivljenost parametrizirane krivulje hitrost, s katero se su£e njena tangenta, medtem ko je torzijska ukrivljenost hitrost sukanja binormale.

Ko govorimo o ukrivljenosti krivulje, imamo v mislih eksijsko ukrivljenost, £e ni drugega dogovora.

Radi bi denirali ²e enostavne sklenjene krivulje, njihovo periodo in orientacijo.

Pri orientaciji bomo najprej stopili korak nazaj in obravnavali ravninske krivulje, podane s parametrizacijo oblike r :I ⊂ R→ R². Bolj splo²no pa lahko deniramo periodi£nost.

Denicija 2.6. Naj bo r : R → Rⁿ gladka parametrizacija krivulje ina naj bo T ∈R. Krivulja, podana s parametrizacijo r, je T-periodi£na, £e velja

r(t+T) =r(t) za vsakt ∈R.

ƒe r ni konstantna preslikava in je T-periodi£na za neki T ̸= 0, potem je krivulja sklenjena.

Denicija 2.7. Enostavna sklenjena krivulja v ravnini je sklenjena ravninska kri- vulja brez samoprese£i²£.

Na ravninske krivulje smo se omejili, ker za njih velja Jordanov izrek o krivu- ljah [12]. Ta izrek, ki je netrivialen rezultat ravninske topologije, nam zagotavlja dobro deniranost orientacije krivulje. Izrek pravi, da ima vsaka enostavna sklenjena krivulja notranjost in zunanjost. Naj bo ta krivulja podana s parametrizacijo r. Bolj natan£no: komplement slike preslikaver je disjunktna unija dveh mnoºic, ki ju ozna£imo z int(r) in ext(r), z naslednjimi lastnostmi:

1. Mnoºica int(r)je omejena, torej vsebovana v krogu z dovolj velikim polmerom.

2. Mnoºica ext(r) je neomejena.

3. Obe mnoºici, int(r) in ext(r), sta povezani, torej za vsaki dve to£ki znotraj iste mnoºice obstaja krivulja, ki je v celoti vsebovana v tej mnoºici, medtem ko vsaka krivulja, ki povezuje to£ko iz mnoºice int(r)s to£ko iz mnoºice ext(r), seka krivuljo r.

Naj bo zdajr naravna parametrizacija ravninske krivulje. Uvedimo oznako t= ṙ. Obstajata dva enotska vektorja, pravokotna na t(s). Deniramo predzna£eno enotsko normalo n_s(s)kot tisti enotski vektor, ki ga dobimo z rotacijo vektorjat(s) za kot π/2 v nasprotni smeri urinega kazalca.

Zaradi dejstva, da ima enostavna sklenjena krivulja notranjost in zunanjost, lahko lo£imo dve orientiranosti krivulje r. Pravimo, da je krivuljar pozitivno orien- tirana, £e njena predzna£ena enotska normala kaºe v notranjost mnoºice int(r). To lahko vedno doseºemo z menjavo parametra t v −t, £e je potrebno. Geometrijsko pozitivna orientacija pomeni, da £e se postavimo na to£ko na krivulji in gledamo v smeri tangentnega vektorja, vidimo int(r) na svoji levi.

Pojem orientacije in enostavne sklenjene ravninske krivulje bomo v naslednjem razdelku posplo²ili na krivulje na ploskvah.

2.2 Parametrizirane ploskve in ukrivljenosti ploskev

Denicija 2.8. Naj bo D odprta mnoºica to£k v ravnini in r : D ⊂ R² → R³ preslikava, ki je razreda C^∞:

r(u, v) = (x(u, v), y(u, v), z(u, v)). (2.1) Ploskev P je mnoºica to£kr(D). Parametrizacija te ploskve je regularna, £e povsod na obmo£juD velja

r_u×r_v ̸= 0,

kjer sta r_u in r_v parcialna odvoda parametrizacije r. Par, ki sestoji iz preslikave r:D→R³ in mnoºice P =r(D), imenujemo regularna parametrizirana ploskev.

Posebnost regularne parametrizirane ploskve je ta, da za vsako to£ko (u₀, v₀) obmo£ja Dobstaja dovolj majhna okolicaD0 te to£ke, ki jo preslikava (2.1) homeo- morfno upodobi na pripadajo£i del ploskveP, kar lahko dokaºemo z uporabo izreka o implicitnih funkcijah.

Od tu naprej predpostavimo, da je na²a parametrizacija ploskve regularna. Naj bo nenotska normala ploskve P, ki jo dolo£a parametrizacija r. Izra£unamo jo kot

n= r_u×r_v

∥r_u×r_v∥.

Normala ploskve v dani to£ki je normala tangentne ravnine v tisti to£ki, ki jo nape- njata vektorjar_u inr_v.

Denirajmo ²e prvo in drugo fundamentalno formo, katerih koeciente upora- bljamo pri ra£unanju ukrivljenosti. Naj bo podana ploskev r(D) in naj bo γ(t) = (u(t), v(t)), t ∈ [0,1], parametrizirana krivulja v D. Prvo fundamentalno formo sre£amo, ko ºelimo izra£unati dolºinoℓ prostorske krivulje

λ(t) =r(γ(t)) =r(u(t), v(t)), ki leºi na ploskvi r(D). Ra£unajmo

ℓ(λ) =

∫︂ 1 0

∥λ̇ (t)∥dt

=

∫︂ 1 0

∥ru(u(t), v(t))·u̇ (t) +rv(u(t), v(t))·v̇ (t)∥dt

=

∫︂ 1 0

(︂√︁(r_uu̇ +r_vv̇ )·(r_uu̇ +r_vv̇ ))︂

(t)dt

=

∫︂ 1 0

(︃√︂

∥r_u∥²u̇² + 2⟨r_u,r_v⟩u̇v̇ +∥r_v∥²v̇² )︃

(t)dt.

ƒe ozna£imo

E =r_u·r_u, F =r_u·r_v, G=r_v·r_v, lahko zapi²emo

ℓ(λ) =

∫︂ 1 0

√︁E(γ(t))u̇ (t)²+ 2F(γ(t))u̇ (t)v̇ (t) +G(γ(t))v̇ (t)²dt.

Izpeljali smo prvo fundamentalno formo ploskve, ki jo z diferenciali zapi²emo kot I =Edu²+ 2F dudv+Gdv²,

E, F in G pa so njeni koecienti. Velja I = ds², to je kvadrat lo£nega elementa krivulje na ploskvi. Pokaºimo, da za regularno parametrizacijo r vedno velja EG− F² >0:

EG−F² = (r_u·r_u)(r_v ·r_v)−(r_u·r_v)²

= (r_u×r_v)·(r_u×r_v)

=∥r_u×r_v∥² >0.

Pri ra£unanju dolºine poti λ(t) = r(u(t), v(t)) smo hkrati izpeljali formulo za tan- gento na to krivuljo:

λ̇ (t0) = r_u(u(t₀), v(t₀))u̇ (t₀) +r_v(u(t₀), v(t₀))v̇ (t₀).

Spomnimo se ²e, da plo²£ino ploskve P izra£unamo kot pl(P) =

∫︂ ∫︂

D

√

EG−F²dudv, kjer jeDdenicijsko obmo£je parametrizacije ploskve. Izraz √

EG−F²dudv =dP je diferencial ploskovne povr²ine.

Preden se lotimo ukrivljenosti ploskve, si oglejmo ukrivljenost krivulj na ploskvi.

Naj bo podana ploskev P in to£ka p na njej. Naj bo r(s) krivulja na ploskvi, parametrizirana z naravnim parametrom, ki gre skozi to£kop. Tangenta na krivuljo je v vsaki to£ki pravokotna na ploskovno normalo n:

n·t= 0.

ƒe opazujemo ploskovno normalo vzdolº krivulje, sta tako n kott odvisna od para- metra s. Zgornjo enakost odvajamo po s in dobimo

ṅ ·t+n·ṫ = 0. (2.2)

Naj bo n_k =ṫ/∥ṫ∥ (pri podpoglavju o krivuljah smo ta vektor ozna£ili z n, a tu je oznaka rezervirana za ploskovno normalo). Upo²tevamo, da je

ṫ = 1 ρn_k, zato je

n·ṫ =n·n_k1 ρ = 1

ρcosϑ.

Pri tem je 1/ρ ukrivljenost krivulje na ploskvi v to£ki p, s ϑ pa smo ozna£ili kot med ploskovno normalo n in n_k, to je glavno normalo krivulje r(s) v to£ki p. To vstavimo v enakost (2.2) in dobimo

cosϑ

ρ =−ṅ ·t.

Ker je ploskovna normala vzdolº krivulje enaka n(s) =n(u(s), v(s)), je ṅ = n_uu̇ + n_vv̇. Sledi

cosϑ

ρ =−(n_uu̇ +n_vv̇ )(r_uu̇ +r_vv̇ ). (2.3) Na tej to£ki denirajmo drugo fundamentalno formo

II=Ldu²+ 2M dudv+N dv², kjer so

L=−n_u·r_u, M =−1

2(n_u·r_v +n_v ·r_u), N =−n_v·r_v. Ko koeciente druge fundamentalne forme vstavimo v enakost (2.3), dobimo

cosϑ

ρ =Lu̇²+ 2M u̇v̇ +N v̇²,

in £e upo²tevamo ²e u̇ = du/ds, v̇ =dv/ds in ds² =Edu²+ 2F dudv+Gdv², lahko zapi²emo

cosϑ

ρ = Ldu²+ 2M dudv+N dv²

Edu²+ 2F dudv+Gdv². (2.4) Na desni strani te enakosti smo dobili ravno kvocient druge in prve fundamentalne forme. Ker ta kvocient vklju£uje koeciente prve in druge fundamentalne forme, bi si ºeleli £im bolj preprosto formulo za njihov izra£un. Natan£neje, radi bi formule za koeciente, kjer ne nastopajo parcialni odvodi normal, temve£ samo parcialni odvodi parametrizacije ploskve. Seveda velja

n·ru = 0,

saj je normala pravokotna na tangentno ravnino. To enakost odvajamo, in vidimo, da je n·r_uu+n_u ·r_u = 0, torej je L = n·r_uu. Iz lastnosti vektorskega produkta sledi, da je∥r_u×r_v∥² =EG−F², zato je

n= r_u×r_v

√EG−F². To dvoje zdruºimo v formulo

L= (r_u×r_v)·r_uu

√EG−F² .

Tako je koecient izraºen samo s parcialnimi odvodi parametrizacije ploskve. Ana- logno vidimo, da je

N = (r_u×r_v)·r_vv

√EG−F² .

Ko odvajamo enakostn·r_u = 0 po v in enakost n·r_v = 0 po u, pa izra£unamo, da je

M = (r_u×r_v)·r_uv

√EG−F² .

ƒe povzamemo:

L=n·ruu, M =n·ruv, N =n·rvv, n= r_u ×r_v

√EG−F².

Oglejmo si podrobneje izrazcosϑ/ρ, ki smo ga izra£unali v (2.4). Ta izraz je odvisen le od razmerja du/dv, to razmerje pa nam dolo£i smer tangente. Iz tu vidimo, da imata krivulji na ploskvi z enako tangento in glavno normalo enako ukrivljenost. Naj bo Σ ravnina, ki vsebuje to£ko p, tangento na²e krivulje in njeno glavno normalo.

Presek ravnineΣin ploskveP je krivulja, ki ima enako ukrivljenost kot na²a za£etna krivulja. Tak presek imenujemo po²evni presek. Presek, ki vsebuje dano to£ko na ploskvi in je pravokoten na tangentno ravnino ploskve v tisti to£ki, pa imenujemo normalni presek. Z normalnim presekom dolo£imo krivuljo na ploskvi, katere glavna normala je vzporedna ploskovni normali. Ker sta obe enotski, velja n = nk ali pa n =−n_k. Kot med njima je torej 0 ali π in posledi£no je cosϑ = ±1. Zdaj lahko deniramo normalno ukrivljenost ploskve.

Denicija 2.9. Naj bo s parametrizacijo r(u, v) = (x(u, v), y(u, v), z(u, v)) podana regularna parametrizirana ploskev P. Izberimo si to£ko p na tej ploskvi in smerni vektore_θ na tangentni ravnini v tej to£ki. Potem je ukrivljenost normalnega preseka ali normalna ukrivljenost ploskve v tej to£ki, glede na izbran smerni vektor, enaka

κ_N(p,e_θ) = 1

R = Ldu²+ 2M dudv+N dv² Edu²+ 2F dudv+Gdv².

Ker je notacija poenostavljena, na prvi pogled ni jasno, kje v formuli se odraºa izbrani vektor e_θ. Skriva se v izbiri krivulje na ploskvi, r(s) = r(u(s), v(s)) in posledi£no v fundamentalnih formah. Kotθpri izbiri smernega vektorja na tangentni ravnini oziroma krivulje na ploskvi prete£e interval [0,2π], pri vsaki izbiri pa lahko izra£unamo pripadajo£o normalno ukrivljenost. Te normalne ukrivljenosti doseºejo svoji ekstremni vrednosti.

Denicija 2.10. Glavni ukrivljenosti ploskve v dani to£ki sta ekstremni vrednosti, ki jih doseºe normalna ukrivljenost, ko izbrani vektor v tangentni ravnini prete£e vse smeri. Ozna£imo ju sκ₁ inκ₂. Vsaki od njiju pripada smerni vektor, pri katerem je ta ukrivljenost doseºena. Ta dva vektorja imenujemo smeri glavnih ukrivljenosti oziroma glavni smeri.

Velja, da sta glavni smeri vedno pravokotni. Za bolj²o predstavo si lahko ravnini glavnih ukrivljenosti, torej ravnini, ki vsebujeta smerna vektorja glavnih ukrivljeno- sti in normalo, ogledamo na sliki 1.

Oglejmo si, kako lahko glavni ukrivljenosti izra£unamo. I²£emo ekstreme funkcije 1

R =Lu̇² + 2M u̇v̇ +N v̇². (2.5) Smer v tangentni ravnini je dolo£ena z u̇ in v̇, zato ju obravnavamo kot iskani ne- znanki. Spomnimo se, da sta v tej formuli odvodau̇ inv̇ glede na naravni parameter s. Zato mora veljati

Eu̇²+ 2F u̇v̇ +Gv̇² = 1.

Slika 1: Ravnini glavnih ukrivljenosti [13].

Vidimo, da smo iskanje ekstremnih vrednosti ukrivljenosti prevedli na vezane ekstreme. V dani to£ki so koecienti prve in druge fundamentalne forme konstantni.

Spomnimo se, da vezane ekstreme i²£emo z Lagrangeevo metodo. Namesto u̇ in v̇ pi²imo xin y ter denirajmo funkcijo

f(x, y) = Lx²+ 2M xy+N y²−λ(Ex²+ 2F xy+Gy²−1).

Parcialna odvoda po x iny morata biti v iskanem ekstremu enaka 0:

2Lx+ 2M y−λ(2Ex+ 2F y) = 0, 2M x+ 2N y−λ(2F x+ 2Gy) = 0.

Pokraj²amo ena£bi z2, prvo ena£bo pomnoºimo zx, drugo zy, se²tejemo in dobimo Lx²+ 2M xy+N y²−λ(Ex²+ 2F xy+Gy²) = 0.

Ker veljaEu̇² + 2F u̇v̇ +Gv̇² = 1, je Lagrangev multiplikator λ enak λ =Lu̇²+ 2M u̇v̇ +N v̇² = 1

R.

Preoblikujmo ena£bi za iskanje stacionarnih to£k pomoºne funkcije f v naslednjo obliko:

(L−λE)x+ (M −λF)y = 0, (M −λF)x+ (N−λG)y = 0.

Vemo, da bo sistem imel neni£elno re²itev, ko bo determinanta matrike det

[︃L−λE M −λF M −λF N −λG ]︃

= 0, (2.6)

saj bo takrat njeno jedro netrivialno. Tako dobimo kvadratno ena£bo za λ:

(EG−F²)λ²−(EN +LG−2F M)λ+ (LN −M²) = 0. (2.7)

Ta ena£ba ima dve realni re²itvi, ki ju ozna£imo s κ₁ = 1

R1 in κ₂ = 1 R2

, in ju imenujemo glavni ukrivljenosti.

Denicija 2.11. Povpre£na ukrivljenost κ je povpre£je normalnih ukrivljenosti κ_N(θ), ko kotθ smernega vektorjae_θ prete£e interval [0,2π]:

κ_H = 1 2π

∫︂ 2π 0

κ_N(θ)dθ. (2.8)

Normalno ukrivljenost lahko izrazimo z glavnima ukrivljenostima, kar nam pove Eulerjev izrek [12].

Izrek 2.12. Naj bo φ kot, ki ga oklepa ravnina Σ normalnega preseka z ravnino normalnega preseka z ukrivljenostjo1/R₁. Potem je ukrivljenost normalnega preseka z ravnino Σ enaka

1 R = 1

R₁ cos²(φ) + 1

R₂ sin²(φ).

To nas pripelje do izraºave povpre£ne ukrivljenosti z glavnima:

κH = κ₁+κ₂

2 . (2.9)

Pokaºimo ²e en na£in za izra£un glavnih ukrivljenosti. Lahko ju izrazimo tudi kot lastni vrednosti matrike ukrivljenosti

B = [︃a b

c d ]︃

, (2.10)

pri £emer moramo a, b, c in d ²e ustrezno dolo£iti. Poi²£imo torej koeciente te matrike. Ena£bo (2.6) lahko zapi²emo kot

det

(︃[︃L M

M N

]︃

−λ

[︃E F

F G

]︃)︃

= 0.

To bo veljalo natanko takrat, ko bo imela ena£ba (︃[︃L M

M N

]︃

−λ

[︃E F

F G

]︃)︃

v= 0 netrivialno re²itev v. To ena£bo preoblikujemo:

[︃L M

M N

]︃

v=λ

[︃E F

F G

]︃

v

[︃E F

F G

]︃−1[︃

L M

M N

]︃

v=λv.

Vidimo, da sta glavni ukrivljenosti lastni vrednosti matrike B =

[︃E F

F G

]︃−1[︃

L M

M N

]︃

= 1

EG−F²

[︃GL−F M GM −F N EM −F L EN −F M ]︃

, (2.11) torej v (2.10) velja

a= GL−F M

EG−F² , b= GM−F N

EG−F² , c= EM −F L

EG−F² , d= EN −F M

EG−F² . (2.12) Preverimo lahko, da z oznakami (2.10) veljaκ_H = ¹₂(a+d) inκ_G =ad−bc, kjer so a, b, c,d podani v (2.12).

Videli bomo, da lahko s pomo£jo matrike ukrivljenosti izra£unamo normalno ukrivljenost v katerikoli smeri v tangentni ravnini. Vpeljimo oznako T_p(M) za tan- gentno ravnino na ploskev M, podano s parametrizacijor(u, v), v to£kip in oznaki F₁ in F₂ za matriki

F1 =

[︃E F

F G

]︃

, F2 =

[︃L M

M N

]︃

. (2.13)

Najprej bi radi videli, da vsak skalarni produkt dveh vektorjev v, w ∈ T_p(M) vklju£uje mnoºenje z matriko F₁, ki vsebuje koeciente prve fundamentalne forme ploskve M. Naj bosta vektorjav inw zapisana v bazi {r_u,r_v}:

v=v₁r_u+v₂r_v, w=w1ru+w2rv. Njun skalarni produkt je enak

v·w=v₁w₁(r_u·r_u) + (v₁w₂ +v₂w₁)(r_u·r_v) +v₂w₂(r_v·r_v)

=v₁w₁E+v₁w₂F +v₂w₁F +v₂w₂G

= [︃v₁

v₂ ]︃T [︃

E F

F G

]︃ [︃

w₁ w₂ ]︃

=vˆ^TF₁wˆ, kjer sta

vˆ = [︃v1

v₂ ]︃

in wˆ = [︃w1

w₂ ]︃

vektorja, ki imata za komponente koeciente vektorjev v in w, ko ju zapi²emo v bazi{r_u,r_v}.Ta ugotovitev je zanimiva, saj takega produkta nismo vajeni iz linearne algebre. Mnoºenje z matriko F₁ se pojavi zaradi dejstva, da na²a baza {r_u,r_v} ni ortonormirana. ƒe bi bila, bi veljalo E = G = 1 in F = 0, zato bi matrika F1

bila kar identiteta in skalarni produkt bi bil enak obi£ajnemu. Zato moramo biti na poljubni ploskvi pazljivi pri uporabi skalarnega produkta in uporabiti izpeljano formulo

⟨v,w⟩=vˆ^TF₁wˆ. (2.14)

Bazo {r_u,r_v} poimenujemo standardna baza, oznaka vˆ pa naj od tu ozna£uje vektor v zapisan v standardni bazi (kot stolpec). Za izpeljavo tega, kako lahko matriko B uporabimo za izra£un normalne ukrivljenosti, potrebujemo naslednje leme.

Lema 2.13. Naj matrikaF₁, denirana z (2.13), vsebuje koeciente prve fundamen- talne forme za ploskev M, podano s parametrizacijo r(u, v). Za dano to£ko p ∈M in vektor v∈T_p(M) velja

vˆ =F₁⁻¹

[︃v·r_u v·r_v ]︃

. Dokaz. Naj bo v=u̇ru+v̇rv. Velja

[︃v·r_u v·rv

]︃

=

[︃u̇r_u·r_u +v̇r_v·r_u u̇ru·rv +v̇rv·rv

]︃

=

[︃u̇E+v̇F u̇F +v̇G

]︃

=

[︃E F

F G

]︃ [︃

u̇ v̇ ]︃

.

Opomba 2.14. Lema 2.13 nam pove, kako lahko s pomo£jo matrike F1 izra£unamo komponente vektorja iz tangentne ravnine, zapisanega v standardni bazi.

Lema 2.15. Naj matrika F₂, podana z (2.13), vsebuje koeciente druge fundamen- talne forme za ploskev M, podano s parametrizacijo r(u, v). Za dano to£ko p ∈M in vektor v∈Tp(M) je normalna ukrivljenost v smeri vektorja v enaka

κN =F2vˆ·vˆ.

Dokaz. Enostavno razpi²emo desno stran enakosti:

F₂vˆ·vˆ =

[︃L M

M N

]︃ [︃

u̇ v̇ ]︃

· [︃u̇

v̇ ]︃

=

[︃Lu̇ +M v̇ M u̇ +N v̇ ]︃

· [︃u̇

v̇ ]︃

=Lu̇²+ 2M u̇v̇ +N v̇²

=κ_N.

ƒe zdruºimo lemi 2.15 in 2.13, vidimo, da je κ_N =F₂F₁⁻¹

[︃v·r_u v·r_v ]︃

·F₁⁻¹

[︃v·r_u v·r_v ]︃

. (2.15)

Normalno ukrivljenost lahko izra£unamo tudi z uporabo matrike ukrivljenosti na slede£i na£in [9]. Za poljuben vektor u=[︁

u₁ u₂]︁T

v tangentni ravnini velja κ_N =u^TF₁⁻¹F₂u=u^TBu. (2.16) Bistvenega pomena pri izpeljavi (2.16) je uporaba formule (2.14).

’e ena pomembna mera ukrivljenosti je Gaussova ukrivljenost.

Denicija 2.16. Gaussova ukrivljenost ploskve v dani to£ki je produkt glavnih ukrivljenosti ploskve v tej to£ki:

κG = 1

R₁R₂. (2.17)

Tako Gaussovo kot povpre£no ukrivljenost lahko z uporabo (2.7) izrazimo s ko- ecienti fundamentalnih form na naslednji na£in:

κ_G = LN −M²

EG−F² , κ_H = 1

2 · EN +GL−2M F

EG−F² . (2.18)

Gaussovo ukrivljenost lahko izra£unamo tudi z ena£bo κ_G= n_u×n_v

r_u×r_v , (2.19)

kjer jer parametrizacija ploskve in nploskovna normala. Dokaz lahko bralec najde v [12, poglavje 8]. Predznak Gaussove ukrivljenosti nam pove nekaj o obliki ploskve:

1. ƒe imata v dani to£ki obe glavni ukrivljenosti enak predznak, torej κ1κ2 >0, potem pravimo, da je ta to£ka elipti£na. Obstaja okolica elipti£ne to£ke, kjer ploskev leºi vsa na eni strani tangentne ravnine.

2. To£ko, v kateri velja κ₁κ₂ < 0, imenujemo hiperboli£na ali sedlasta to£ka. V taki to£ki ima ploskev obliko sedla.

3. Naj bo v dani to£ki κ₁κ₂ = 0. ƒe sta obe glavni ukrivljenosti enaki 0, je to£ka planarna. ƒe je ena izmed glavnih ukrivljenosti neni£elna, pravimo, da je to£ka paraboli£na.

Oglejmo si ²e geodetsko ukrivljenost krivulje na dani ploskvi P. Naj bor njena naravna parametrizacija. Vemo ºe, da je ṙ enotska vektorska funkcija, ki je seveda tangentna na ploskev in pravokotna na ploskovno normalo n ploskve P. Torej so vektorske funkcije n, ṙ in n× ṙ paroma pravokotne. Spet, ker je r naravna parametrizacija, je tudi vektorska funkcijar¨pravokotna naṙ, zato jo lahko zapi²emo kot linearno kombinacijo n inn×ṙ,

r¨ =κNn+κgn×ṙ, (2.20) kjer je κ_N ºe poznana normalna ukrivljenost.

Denicija 2.17. Skalar κ_g iz enakosti (2.20) je geodetska ukrivljenost krivulje, podane s parametrizacijor.

Izpeljano si lahko ogledamo na sliki 2. Na kratko si oglejmo ²e geodetske krivulje.

Denicija 2.18. Krivulja na ploskvi P, podana s parametrizacijo r, je geodetska krivulja ali geodetka, £e je za vsak parameter t vektor ¨(t)r ni£eln ali pa pravokoten na tangentno ravnino ploskve v to£ki r(t), torej vzporeden enotski normali.

Slika 2: Geodetska ukrivljenost krivulje [12].

Lahko si predstavljamo, da smo to£ka na ploskvi. Geodetka bi bila za nas tista krivulja, ki bi jo dojeli kot ravno [12]. Najkraj²a pot med dvema to£kama na ploskvi je na primer vedno geodetka. ’e ena primerjava je voºnja po ravni cesti mi jo dojemamo kot ravno, £eprav povr²je zemlje ni ravno. ƒe je cesta podana s parametrizacijo r, njen pospe²ek r¨ni ni£eln, vendar jo dojemamo kot ravno, ker je tangentna komponenta r¨ ni£elna (r¨ je pravokoten na sfero, ki predstavlja povr²je zemlje).

Geodetka ima vedno konstantno hitrost in njena naravna reparametrizacija bo

²e vedno geodetka, a ni res, da bo vsaka reparametrizacija geodetke ²e vedno geo- detka [12]. Brez dokaza navedimo, da je naravno parametrizirana krivulja geodetka natanko takrat, ko je njena geodetska ukrivljenost konstantno enaka ni£ (bralec lahko dokaz najde v [12, poglavje 9]).

Trditev 2.19. Vsak izsek premice je geodetka.

Dokaz. Vsako premico lahko parametriziramo kot r(t) =a+bt,

kjer sta ain b konstantna vektorja. O£itno veljar¨(t) = 0.

Spomnimo se ²e orientabilnih in orientiranih ploskev.

Denicija 2.20. Ploskev P ⊂R³ je orientabilna, £e obstaja zvezno vektorsko polje n naR³, denirano na P, tako da za vsako to£ko p ∈P obstaja enotski vektor n_p pravokoten na tangentno ravnino na ploskev P v to£kip.

ƒe je ploskev orientabilna, imamo dve izbiri, n in −n. Oglejmo si ²e denicijo enostavne sklenjene krivulje na ploskvi in njene orientacije.

Denicija 2.21. Krivulja γ(t) = r(u(t), v(t)) na ploskvi r : U →R³ je enostavna sklenjena krivulja s periodo T, £e jeπ(t) = (u(t), v(t))enostavna sklenjena krivulja v R² s periodo T in velja int(π) ⊂ U. Krivulja γ je pozitivno orientirana, £e je π pozitivno orientirana. Deniramo ²e

int(γ) = r(int(π)).

3 Ra£unanje ukrivljenosti s pomo£jo aproksimacije z Bézierjevimi ploskvami

Eden izmed na£inov ra£unanja ukrivljenosti prostorske triangulacije je aproksima- cija triangulacije z dovolj gladko ploskvijo, na kateri lahko uporabimo standardne formule za ra£unanje poljubno izbrane ukrivljenosti. Klju£nega pomena je seveda dobra aproksimacija dane triangulacije. V naslednjem podpoglavju je predstavljena metoda aproksimacije z Bézierjevimi ploskvami iz tenzorskega produkta. Preden za£nemo z osnovami Bézierjevih krivulj in ploskev denirajmo prostorsko triangula- cijo.

Denicija 3.1. Prostorska triangulacijaKnad to£kami{xi}je mnoºica trikotnikov v R³, za katero velja

ˆ Ogli²£a trikotnikov so to£ke xi.

ˆ Notranjosti poljubnih dveh trikotnikov se ne sekata.

ˆ ƒe imata dva trikotnika neprazen presek, imata skupno ogli²£e ali pa stranico.

3.1 Uvod v Bézierjeve krivulje in Bézierjeve ploskve

Bézierjeve krivulje in ploskve so zaradi svojih lastnosti zelo uporabne na podro-

£ju ra£unalni²ko podprtega geometrijskega oblikovanja. Pri njihovi konstrukciji so klju£ni Bernsteinovi bazni polinomi, saj so te krivulje in ploskve denirane kot nji- hove linearne kombinacije.

Denicija 3.2. Za i= 0,1, . . . nje i-ti Bernsteinov bazni polinom stopnje n deni- ran kot

B_iⁿ(t) = (︃n

i )︃

tⁱ(1−t)ⁿ⁻ⁱ.

Primer 3.3. Oglejmo si nekaj primerov Bernsteinovih baznih polinomov:

1. n= 0: B₀⁰(t) = 1,

2. n= 1: B₀¹(t) = 1−t, B₁¹(t) =t,

3. n= 2: B₀²(t) = (1−t)², B₁²(t) = 2t(1−t),B₂²(t) =t².

S kraj²im ra£unom lahko preverimo, da zado²£ajo rekurzivni zvezi B_iⁿ(t) = (1−t)B_iⁿ⁻¹(t) +tBⁿ⁻¹_i−1(t),

in da za odvod velja

(B_iⁿ(t))^′ =n(Bⁿ⁻¹_i−1(t)−B_iⁿ⁻¹(t)).

Na²tejmo ²e nekaj lastnosti Bernsteinovih baznih polinomov:

ˆ Bⁿ_i(0) =B_iⁿ(1) = 0za vsak i= 1,2, . . . , n−1in n≥2,

ˆ Bⁿ₀(0) = 1,B₀ⁿ(1) = 0 za vsakn ≥1,

ˆ Bⁿ_n(0) = 0,B_nⁿ(1) = 1 za vsakn ≥1.

Ozna£imo vektorski prostor polinomov stopnje ≤n s Pn=Lin{1, t, t², . . . , tⁿ}. Trditev 3.4. Bernsteinovi bazni polinomi tvorijo bazo za Pn.

Dokaz. Za dokaz trditve uporabimo indukcijo.

1. Preverimo, da trditev velja za n = 0. Res, B⁰₀(t) = 1. Pokaºimo ²e, da velja tudi za n = 1. Bernsteinova polinoma sta enaka B₀¹(t) = 1−t in B₁¹(t) = t. Tedaj lahko zapi²emo

[︃1 t ]︃

= [︃1 1

0 1 ]︃ [︃

B₀¹(t) B₁¹(t) ]︃

, kar dokazuje veljavnost trditve za n= 1.

2. Recimo, da izrek velja za stopnjo n −1. Dovolj je dokazati, da so B_iⁿ, i = 0,1, . . . , n, linearno neodvisni. Ker jih je n + 1, bo to pomenilo, da tvorijo bazo. šelimo dokazati, da za vsak t velja

n

∑︂

i=0

α_iB_iⁿ(t) = 0 ⇐⇒ α₀ =α₁ =. . .=α_n = 0. (3.1)

ƒe velja αi = 0 za vsak i, je o£itno tudi vsota ∑︁n

i=0αiB_iⁿ(t) enaka 0. Naj bo

n

∑︂

i=0

α_iB_iⁿ(t) = 0. (3.2)

V to enakost najprej vstavimo t = 0. Z uporabo lastnosti Bernsteinovih polinomov in upo²tevanjem, da je n ≥ 2, dobimo α₀ = 0. Vstavimo ²e t = 1 in vidimo, da je prav tako α_n = 0. Odvajamo enakost (3.2) in dobimo

(︄ _n

∑︂

i=0

α_iB_iⁿ(t) )︄′

= 0.

Odvod na levi strani enakosti se poenostavi v

n

∑︂

i=0

α_i(B_iⁿ(t))^′ =

n

∑︂

i=0

α_in(B_i−1ⁿ⁻¹(t)−B_iⁿ⁻¹(t))

=n

n

∑︂

i=0

α_iB_i−1ⁿ⁻¹(t)−n

n

∑︂

i=0

α_iB_iⁿ⁻¹(t)

=n

n−1

∑︂

i=0

(α_i+1−α_i)B_iⁿ⁻¹(t),

kjer smo v prvi vsoti zamaknili indeks in upo²tevaliB₋₁ⁿ⁻¹(t) = 0terB_nⁿ⁻¹(t) = 0. Po indukcijski predpostavki so B_iⁿ⁻¹, i= 0, . . . , n−1, linearno neodvisni, zato mora veljatiαi+1−αi = 0 za vsak i= 0,1, . . . , n−1, od koder z upo²te- vanjemα₀ = 0 dobimo

α0 =α1 =α2 =. . .=αn = 0, kar dokazuje trditev.

Iz dokaza je razvidno, kako odvajamo linearno kombinacijo Bernsteinovih baznih polinomov.

Posledica 3.5. Velja naslednja formula:

(︄ _n

∑︂

i=0

α_iB_iⁿ(t) )︄′

=n

n−1

∑︂

i=0

(α_i+1−α_i)B_iⁿ⁻¹(t).

Denicija 3.6. Naj bodo dane to£ke b₀, b₁, . . . , b_n v R^d. Bézierjeva krivulja je parametri£no podana krivuljapn: [0,1]→R^d, denirana s predpisom

pn(t) =

n

∑︂

i=0

biB_iⁿ(t).

To£ke b₀, b₁, . . . , b_n imenujemo kontrolne to£ke, poligon, ki jih povezuje, pa kontrolni poligon. Iz lastnosti Bernsteinovih baznih polinomov

n

∑︂

i=0

B_iⁿ(t) = 1, B_iⁿ(t)≥0, za vsakt ∈[0,1],

sledi, da je za vsak t to£ka na Bézierjevi krivulji konveksna kombinacija kontrolnih to£k, zato krivulja leºi v njihovi konveksni ovojnici. Analogno lahko deniramo Bézierjevo ploskev.

Denicija 3.7. Bézierjeva ploskev (iz tenzorskega produkta) je parametri£no po- dana polinomska ploskev

p_m,n : [0,1]×[0,1]→R³, p_m,n(u, v) =

m

∑︂

i=0 n

∑︂

j=0

b_i,jB_i^m(u)B_jⁿ(v), kjer sob_i,j kontrolne to£ke, ki tvorijo kontrolno mreºo.

Slika 3: Primer Bézierjeve ploskve iz tenzorskega produkta zam =n = 3 [6].

Ploskev p_m,n je sled Bézierjeve krivulje stopnje m, ki se giblje v prostoru in spreminja svojo obliko. V vsakem trenutku jo dolo£a mnoºicam+ 1kontrolnih to£k, od katerih se vsaka giblje po svoji Bézierjevi krivulji stopnje n. Primer Bézierjeve ploskve za m=n= 3 si lahko ogledamo na sliki 3.

Dve lastnosti, zaradi katerih so Bézierjeve ploskve tako primerne za aproksimacijo diskretizirane ploskve, sta naslednji:

ˆ Interpolacija to£kb_0,0, b_0,n,b_m,0 in b_m,n.

ˆ Bézierjeva ploskev leºi v konveksni ovojnici svojih kontrolnih to£k.

To nam zagotavlja lokalno kontrolo nad obliko ploskve.

Za izra£un to£k tako na Bézierjevih krivuljah kot na ploskvah pri danem parame- tru uporabljamo de Casteljaujev algoritem. ƒe za Bézierjevo ploskev velja m = n, se algoritem poenostavi v spodnjo obliko:

Input: b_i,j i, j = 0,1, . . . , n b^0,0_i,j =bi,j i, j = 0,1, . . . , n for r = 1 :n do

for i, j = 0,1, . . . , n−r do b^r,r_i,j(u, v) =[︁

1−u u]︁

[︃b^r−1,r−1_i,j b^r−1,r−1_i,j+1 b^r−1,r−1_i+1,j b^r−1,r−1_i+1,j+1

]︃ [︃

1−v v

]︃

end for end for

return b^n,n_0,0(u, v)

Opomba 3.8. Pri mnoºenju matrike in vektorjev v algoritmu gre za mnoºenje po komponentah.

Na sliki 4 vidimo gra£no ponazoritev de Casteljaujevega algoritma za m = n. Preverimo, da s tem algoritmom res dobimo ustrezno to£ko na ploskvi.

Trditev 3.9. Velja enakost b^r,r_i,j(u, v) =

r

∑︁

k=0 r

∑︁

ℓ=0

b_i+k,j+ℓB_k^r(u)B_ℓ^r(v).

Dokaz. Trditev dokaºemo z indukcijo na r, pri £emer upo²tevamo algoritem in re- kurzivno zvezo za Bernsteinove polinome.

Slika 4: Prikaz de Casteljaujevega algoritma za ploskve iz tenzorskega produkta [6].

1. Preverimo, da trditev velja za r= 0:

0

∑︂

k=0 0

∑︂

ℓ=0

b_i+k,j+ℓB_k⁰(u)B_ℓ⁰(v) = b_i,jB₀⁰(u)B₀⁰(v)

=b_i,j =b^0,0_i,j.

2. Predpostavimo, da trditev velja zar−1, in dokaºimo, da velja zar. Z uporabo algoritma dobimo

b^r,r_i,j(u, v) = b^r−1,r−1_i,j (1−u)(1−v) +b^r−1,r−1_i,j+1 (1−u)v +b^r−1,r−1_i+1,j u(1−v) +b^r−1,r−1_i+1,j+1uv.

Z upo²tevanjem indukcijske predpostavke vidimo, da lahko to izrazimo z

r−1

∑︂

k=0 r−1

∑︂

ℓ=0

b_i+k,j+ℓB_k^r−1(u)B_ℓ^r−1(v)(1−u)(1−v)

+

r−1

∑︂

k=0 r−1

∑︂

ℓ=0

b_i+k,j+ℓ+1B^r−1_k (u)B^r−1_ℓ (v)(1−u)v

+

r−1

∑︂

k=0 r−1

∑︂

ℓ=0

b_i+k+1,j+ℓB^r−1_k (u)B^r−1_ℓ (v)u(1−v)

+

r−1

∑︂

k=0 r−1

∑︂

ℓ=0

bi+k+1,j+ℓ+1B_k^r−1(u)B_ℓ^r−1(v)uv.

ƒe izpostavimo ustrezne £lene, uvedemo nova indeksa in upo²tevamo rekur- zivno zvezo za Bernsteinove bazne polinome, lahko zgornji izraz preoblikujemo

v

r

∑︂

k=0 r

∑︂

ℓ=0

b_i+k,j+ℓB_k^r−1(u)(1−u)(B_ℓ^r−1(v)(1−v) +B_ℓ−1^r−1(v)v) +

r

∑︂

k=0 r

∑︂

ℓ=0

b_i+k,j+ℓB^r−1_k−1(u)u(B_ℓ^r−1(v)(1−v) +B^r−1_ℓ−1(v)v)

=

r

∑︂

k=0 r

∑︂

ℓ=0

b_i+k,j+ℓB^r_ℓ(v)(B_k^r−1(u)(1−u) +B_k−1^r−1(u)u)

=

r

∑︂

k=0 r

∑︂

ℓ=0

bi+k,j+ℓB^r_ℓ(v)B_k^r(u).

Sedaj v enakost iz trditve 3.9 vstavimo r = n in i = j = 0 in neposredno sledi pravilnost delovanja algoritma:

b^n,n_0,0(u, v) =

n

∑︂

k=0 n

∑︂

ℓ=0

b_k,ℓB_kⁿ(u)B_ℓⁿ(v) =p_n,n(u, v).

Oglejmo si ²e pravilo za odvajanje Bézierjevih ploskev. Velja

∂^r+sp_m,n

∂u^r∂v^s (u, v) = m!n!

(m−r)!(n−s)!

m−r

∑︂

i=0 n−s

∑︂

j=0

∆^r,sb_i,jB^m−r_i (u)B_j^n−s(v).

Pri tem je operator ∆deniran na naslednji na£in:

∆^1,0b_i,j =b_i+1,j−b_i,j,

∆^0,1b_i,j =b_i,j+1−b_i,j,

∆^r,sb_i,j = ∆^r−1,sb_i+1,j−∆^r−1,sb_i,j

= ∆^r,s−1b_i,j+1−∆^r,s−1b_i,j.

3.2 Aproksimacija prostorske triangulacije z Bézierjevo plo- skvijo po metodi najmanj²ih kvadratov

Dano prostorsko triangulacijo ºelimo lokalno aproksimirati z Bézierjevo ploskvijo.

Najprej si moramo izbrati stopnji m in n Bézierjeve ploskve. Ker pri ra£unanju ukrivljenosti nastopajo odvodi do drugega reda, si izberemo m = n = 2, saj na ta na£in obdrºimo dovolj informacij o triangulaciji, ra£unanje pa ostane enostavno.

Iskana ploskev bo torej oblike p2,2(u, v) =

2

∑︂

i=0 2

∑︂

j=0

bi,jB_i²(u)B_j²(v).

Dolo£iti moramo kontrolne to£keb_i,j, ki jih pi²emo kot vrstice. Denimo, da triangu- lacijo dolo£ajo to£ke t_k, k = 0,1, . . . , K. Predpostavimo, da vsaki to£ki t_k pripada

par parametrov(u_k, v_k). šelimo si, da bi veljalot_k =p_2,2(u_k, v_k), kar lahko zapi²emo kot

tk =[︁

B₀²(u_k)B₀²(v_k), . . . , B₂²(u_k)B₂²(v_k)]︁

⎡

⎢

⎣ b_0,0

...

b_2,2

⎤

⎥

⎦,

⎡

⎢

⎣ b_0,0

...

b_2,2

⎤

⎥

⎦∈R^9×3.

Zaradi preglednosti je ta zapis okraj²an in iz njega ni razviden vrstni red kom- ponent vektorjev. Razporedimo jih lahko tako, da paru (i, j) priredimo indeks ℓ= (m+ 1)·i+j+ 1 = 3i+j + 1, i, j = 0,1,2.

ƒe zdruºimo vseh K+ 1 ena£b, dobimo sistem linearnih ena£b

⎡

⎢

⎢

⎢

⎢

⎢

⎢

⎣ t₀

...

...

t_K

⎤

⎥

⎥

⎥

⎥

⎥

⎥

⎦

=

⎡

⎢

⎢

⎢

⎢

⎢

⎣

B₀²(u₀)B₀²(v₀) . . . B₂²(u₀)B₂²(v₀)

... ... ...

B₀²(u_K)B₀²(v_K) . . . B₂²(u_K)B₂²(v_K)

⎤

⎥

⎥

⎥

⎥

⎥

⎦

⎡

⎢

⎢

⎢

⎣ b_0,0

...

b_2,2

⎤

⎥

⎥

⎥

⎦ ,

kjer imamo ve£ ena£b kot neznank. Tak sistem imenujemo predolo£en sistem in ga v matri£ni obliki zapi²emo kot

T =M B, (3.3)

kjer sta matrikiT ∈R^(K+1)×3 in B ∈R^9×3 enaki

T = [T₁, T₂, T₃] =

⎡

⎢

⎢

⎢

⎢

⎢

⎢

⎣ t₀

...

...

t_K

⎤

⎥

⎥

⎥

⎥

⎥

⎥

⎦

, B = [B₁, B₂, B₃] =

⎡

⎢

⎢

⎢

⎢

⎢

⎢

⎢

⎢

⎢

⎢

⎢

⎢

⎣ b_0,0 b_0,1 b_0,2 b_1,0 b_1,1 b_1,2 b_2,0 b_2,1 b_2,2

⎤

⎥

⎥

⎥

⎥

⎥

⎥

⎥

⎥

⎥

⎥

⎥

⎥

⎦ ,

M ∈R^(K+1)×9 pa je matrika

M =

⎡

⎢

⎢

⎢

⎢

⎢

⎣

B₀²(u0)B₀²(v0) . . . B₂²(u0)B₂²(v0)

... ... ...

B₀²(u_K)B₀²(v_K) . . . B₂²(u_K)B₂²(v_K)

⎤

⎥

⎥

⎥

⎥

⎥

⎦ .

Za vsakℓ= 1,2,3moramo torej re²iti predolo£en sistemTℓ =M Bℓ sK+1ena£bami in9 neznankami, katerega formalna re²itev po metodi najmanj²ih kvadratov je

B_ℓ= (M^TM)⁻¹M^TT_ℓ, ℓ = 1,2,3, oziroma B = (M^TM)⁻¹M T.

Slika 5: Prikaz vektorjevM B_ℓ in T_ℓ glede na re²itev normalnega sistema.

Spomnimo se, kako pridemo do te re²itve, in si oglejmo nekaj moºnosti za njen numeri£ni izra£un. Fiksirajmoℓ ∈ {1,2,3}.Re²itev po metodi najmanj²ih kvadratov je tisti vektor B_ℓ, ki minimizira razdaljo ∥M B_ℓ −T_ℓ∥₂. Vemo, da vsi vektorji M B_ℓ leºijo v podprostoru Im(M). Najbliºji element v Im(M)elementu T_ℓ v drugi normi je pravokotna projekcija vektorja T_ℓ na Im(M), kot lahko vidimo na sliki 5. I²£emo torej vektor B_ℓ, za katerega bo veljalo M B_ℓ−T_ℓ ⊥ Im(M). ƒe so stolpci matrike M baza za Im(M), dobimo pogoj M^T(M B−T) = 0, oziroma

M^TM B_ℓ =M^TT_ℓ, (3.4)

kar imenujemo normalni sistem. Re²itev lahko numeri£no izra£unamo na ve£ na£i- nov, ena moºnost je uporaba QR razcepa, kjer matriko M zapi²emo kot produkt matrik QR, kjer jeQmatrika velikost(K+ 1)×9z ortonormirani stolpci inR zgor- nje trikotna matrika s pozitivnimi diagonalnimi elementi. Potem lahko normalni sis- tem (3.4) zapi²emo kot (QR)^T(QR)B_ℓ = (QR)^T_ℓ, kar preoblikujemo v RB_ℓ =Q^TT_ℓ. Lahko bi se zgodilo, da matrika M ne bi bila maksimalnega ranga9. V tem primeru lahko uporabimo t. i. psevdoinverz matrike M, ki ga ozna£imo z M⁺.

Glavni problem, ki nastopi, pa je, da v praksi nimamo podanih prej omenjenih parametrov (u_k, v_k). Opi²imo enega izmed preprostej²ih na£inov, kako jih konstru- iramo. Najprej to£ke t_k = (x_k, y_k, z_k) projiciramo na ravnino, tako da izpustimo njihovo z-koordinato. Dobimo par (x_k, y_k). Nato koordinate transformiramo, tako da leºijo na intervalu [0,1]×[0,1], na koncu ²e dolo£imo u_k = x_k in v_k = y_k. Ta na£in sicer ni optimalna izbira parametrov, ²e posebej ne v primeru, ko imamo dele domene v obmo£ju [0,1]×[0,1], v katerih ni to£k, kot je prikazano na sliki 6. Za bolj²o parametrizacijo si ºelimo £im manj²i pravokotnik, v katerem leºijo to£ke. Ta problem lahko re²imo tako, da najdemo osi najmanj²e elipse, ki vsebuje vse to£ke, in za domeno izberemo najmanj²i pravokotnik, ki to elipso vsebuje. Druga moºna re²itev bi bila, da bi poiskali kak²no drugo smer projekcije, ki ni xy-ravnina.

Re²itev po metodi najmanj²ih kvadratov nam da v praksi dobro aproksimacijo vhodnih podatkov. Teºava, ki se pri tem pojavi, pa je, da so podane prostorske triangulacije v praksi pridobljene s pomo£jo digitalnih naprav, pri katerih veliko- krat pride do ²umov v podatkih. Ta re²itev seveda sledi tudi tem ²umom, £emur se ºelimo izogniti. Eden izmed pokazateljev prisotnosti ²uma je pretirana zvitost kontrolne mreºe, zato postavimo ²e dodatne omejitve, ki omejujejo t.i. zasuk (angl.

Slika 6: Primer porazdelitve projiciranih to£k triangulacije.

twist) kontrolne mreºe. Tega merimo z me²anim odvodom parametrizacije ploskve

∂²

∂u∂v. Zasuk bikvadrati£ne Bézierjeve ploskve je Bézierjeva ploskev stopnje (1,1) s koecienti∆^1,1b_i,j, kjer∆^1,1 meri deviacijo vsakega ²tirikotnika kontrolne mreºe od paralelograma. Res, naj bo q_i,j £etrta to£ka paralelograma, ki ga dolo£ajo to£ke b_i,j, b_i+1,j, b_i,j+1. Paralelogram je dolo£en s pogojem

q_i,j−b_i+1,j =b_i,j+1−b_i,j.

ƒe uporabimo zapis∆^1,1b_i,j, dobimo, da je

∆^1,1b_i,j =b_i+1,j+1−q_i,j.

Vidimo, da operator ∆^1,1 na to£ki b_i,j meri deviacijo vsakega ²tirikotnika v mreºi kontrolnih to£k od paralelograma. Zato prej deniranemu predolo£enemu sistemu (3.3) dodamo ena£be ∆^1,1b_i,j = 0, i, j = 0,1, kar se prepi²e v

b_i+1,j+1−b_i,j+1−b_i+1,j+b_i,j = 0, i, j = 0,1. (3.5) Za vsako komponento neznanih kontrolnih to£k dobimo dodatne ²tiri ena£be, zato vektor T_ℓ raz²irimo v

T˜︁_ℓ =

⎡

⎢

⎢

⎢

⎢

⎣ T_ℓ

0 0 0 0

⎤

⎥

⎥

⎥

⎥

⎦

, ℓ= 1,2,3.

Vektor B_ℓ seveda ostaja enak, spremeni pa se matrika predolo£enega sistema, in sicer se raz²iri v matriko

M˜︂= [︃M

S ]︃

∈R^(K+5)×9, kjer je S ∈R^4×9 enaka

S =

⎡

⎢

⎢

⎣

1 −1 0 −1 1 0 0 0 0

0 1 −1 0 −1 1 0 0 0

0 0 0 1 −1 0 −1 1 0

0 0 0 0 1 −1 0 −1 1

⎤

⎥

⎥

⎦ .