Własny silnik graficzny. Część II: detekcja kolizji.

30.11.2010 - Robert Kraus
TrudnośćTrudnośćTrudnośćTrudność

Badanie przecięcia promień-sfera

1
2
3
4
5
6
7
8
9
10
11
inline float intersect(sphere::i S, ray::i r) {
  vec3 v = r.origin - S.center;
  float b = dotProd(v, r.aheadDir);
  float d2 = b*b - dotProd(v,v) + S.radius*S.radius;
  if (d2 < 0.0f) return undefined;
  float d = sqrt(d2);
  float t1 = -b-d, t2 = -b+d;
  if (t1 > epsilon5) return t1;
  if (t2 > epsilon5) return t2;
  return undefined;
}

Równanie sfery:

$$(x-c_x)^2 + (y-c_y)^2 + (z-c_z)^2 = r^2$$
Wstawiając $ origin + aheadDir * t $ za $ (x,y,z) $ otrzymujemy równanie:
$$((o_x + d_x*t) -c_x)^2 + ((o_y + d_y*t)-c_y)^2 + ((o_z + d_z*t)-c_z)^2 = r^2$$
Równanie to można zapisać (grupując wyrazy ze względu na $ t $) jako: $ At^2 + Bt + C = 0 $, gdzie:
$$A = d_x^2 + d_y^2 + d_z^2$$
$$B = 2(d_x (o_x - c_x) + d_y (o_y - c_y) + d_z (o_z - c_z))$$
$$C = (o_x-c_x)^2 + (o_y-c_y)^2 + (o_z-c_z)^2 - r^2$$

Jeśli $ aheadDir $ ma długość $ 1 $ (a dbamy o to, żeby tak było), to $ A = 1 $. Przyspiesza to nieco obliczenia, a pierwiaski równania $ At^2 + Bt + C = 0 $ możemy obliczyć ze wzorów:

$$ t_1 = \frac{-B - \sqrt{B^2 - 4C}}{2} $$
$$ t_2 = \frac{-B + \sqrt{B^2 - 4C}}{2} $$
W wierszach nr 2-4 obliczamy $ d2 = B^2 - 4C $.
W wierszu nr 5 sprawdzamy czy $ d2 < 0 $, jeśli tak jest, to równanie $ At^2 + Bt + C = 0 $ nie ma rozwiązań. W przeciwnym wypadku przechodzimy do wierszy nr 6-7 i obliczamy $ d = \sqrt{d2} = \sqrt{B^2 -4C} $, a nastepnie $ t1 $ i $ t2 $.

Szukamy najbliższego przecięcia, a więc interesuje nas ta wartość $ t $, która jest większa od zera (jeśli $ t $ jest mniejsza od zera, to punkt przecięcia jest za źródłem, a nie przed źródłem promienia względem kierunku ruchu). Istnieją teraz trzy możliwe sytuacje rozpatrywane w wierszach 8-10.

1) $ t1 $ oraz $ t2 $ są większe od zera, a bliższe przecięcie opisuje $ t1 $.

2) Tylko $ t2 $ jest większe od zera, a więc opisuje interesujące nas przecięcie $ t1 $.

3) $ t1 $ oraz $ t2 $ są mniejsza od zera, a więc brak przecięć.

Badanie przecięcia promień-trójkąt

Do badania przecięcia promienia z trójkątem będziemy potrzebowali funkcji obliczającej pole trójkąta reprezentowanego przez trzy punkty. Uzyjemy do tego implementacji wzoru Herona.

1
2
3
4
5
6
7
inline float area(vec3::i v1, vec3::i v2, vec3::i v3) {
  float a = length(v2 - v1);
  float b = length(v3 - v1);
  float c = length(v3 - v2);
  float p = (a + b + c) / 2.0f;
  return sqrt(p*(p-a)*(p-b)*(p-c));
}

Badanie przecięcia odbywa się w dwóch etapach. Pierwszy etap polega na znalezieniu punktu przecięcia promienia z płaszczyzną, na której leży trójkąt.

Przypomnijmy, że
$ dotProd(N,(x,y,z)) + d = 0 $ to równanie płaszczyzny trójkąta
$ N $ to wektor normalny do powierzchni trójkąta
$ (x,y,z) $ to dowolny punkt przestrzeni
$ origin + aheadDir \cdot a $ oznacza punkty na prostej po której podróżuje promień

Podstawiamy $ origin + aheadDir \cdot a $ za $ (x,y,z) $ otrzymując równanie

$$dotProd(N, origin + aheadDir \cdot t) + d = 0,$$
które rozwiązujemy ze względu na $ t $. Z rozwiązania otrzymujemy, że
$$t = \frac{-d - dotProd(N, origin)}{dotProd(N, aheadDir)}.$$
Jeśli $ t $ jest większe od zera oraz nie jest nieskończonością, to nasz promień trafia w płaszczyznę, na której leży trójkąt. Poniżej kod procedury badającej przecięcia, etap pierwszy wykonywany jest w wierszach 2-5 procedury $ intersect $.

1
2
3
4
5
6
7
8
9
10
11
12
13
inline float intersect(triangle::i T, ray::i r) {
  float t = (-T.d - dotProd(T.N, r.origin)) / dotProd(T.N, r.aheadDir);
  if ( !(t > epsilon5 && t < infinity) )
    return undefined;
  vec3 ip = r.origin + r.aheadDir * t;
  float abc = area(T.v1 , T.v2 , T.v3);
  float pbc = area( ip  , T.v2 , T.v3);
  float apc = area(T.v1 ,  ip  , T.v3);
  float abp = area(T.v1 , T.v2 ,  ip );
  if ( pbc + apc + abp > abc + epsilon3 )
    return undefined;
  return t;
}

Mając punkt przecięcia promienia z płaszczyzną trójkąta (jest to punkt $ ip $ obliczony w wierszu nr 5) możemy przystąpić do drugiego etapu badania przeciącia promienia z trójkętem. Drugi etap polega na sprawdzeniu czy punkt przecięcia promienia z płaszczyzną leży wewnątrz trójkąta.

Niech $ A $, $ B $, $ C $ oznaczają wierzchołki $ v1 $, $ v2 $, $ v3 $ naszego trójkąta, a $ P $ oznacza punkt przecięcia $ ip $. Aby sprawdzić czy $ ip $ leży wewnątrz trójkąta należy sprawdzić czy suma pól trójkątów $ PBC $, $ APB $, $ ABP $ jest równa polu trójkąta $ ABC $ (wtedy punkt $ ip $ leży wewnątrz trójkąta), czy jest większa od pola trójkąta $ ABC $ (wtedy punkt $ ip $ leży poza trójkątem). Tutaj też jednak narażeni jesteśmy na pułapkę. Z powodu niedokładnej reprezentacji liczb rzeczywistych nie jesteśmy w stanie sprawdzać równości wyżej wymienionych pól trójkątów, więc sprawdzamy czy suma pól trójkątów $ PBC $, $ APB $ oraz $ ABP $ jest większa od pola trójkąta $ ABC $ zwiększonego o pewną małą wartość. Etap drugi wykonywany jest w wierszach 6-12 procedury $ intersect $.

Błędy numeryczne

Gdy promień trafi w pewien punkt $ ip $, a następnie wygenerujemy np. promień odbity, to źródłem tego promienia będzie właśnie punkt $ ip $. Procedura dokonująca śledzenia promienia dokona testu przecięcia również z obiektem, od którego ten promień się odbija. W czym problem? Odległość punktu $ ip $ od nowego punktu przecięcia powinna wynosić zero, niestety z powodu niedokładnej reprezentacji liczb rzeczywistych w komputerze na ogół tak nie będzie. Będzie to wartość bardzo bliska zeru, ale jednak nie zero. Procedura śledząca promień uzna więc, że obiekt, od którego promień się odbija jest najbliższym obiektem przecinającym nasz odbity promień, co oczywiście jest błędem, który prowadzi do sytuacji, w której z obliczeń wynika, że powierzchnia sama sobie zasłania źródło światła mimo, że jest to nieprawdą. Skutkuje to błędnym niedoświetleniem różnych punktów sceny objawiajacego się ziarnistym szumem. Oto efekt omawianych błędów obliczeniowych:

Jak sobie z tym poradzić?

Śledzony promień (niebieski) trafił w punkt $ ip $, po czym wygenerowaliśmy promień odbity (czerwony) i rozpoczynamy dla niego procedurę śledzenia. Jeśli w wyniku śledzenia promienia odbitego otrzymamy nowy punkt przecięcia, który znajduje się w obszarze o promieniu $ err $ (rozmiar błędu wynikającego z niedokładnej reprezentacji liczb rzeczywistych) na powyższej ilustracji, to najpewniej jest to wynik błędów numerycznych. Stąd akceptujemy tylko punkty leżące poza obszarem o promieniu $ eps $ odrzucając wszystkie błędne punkty przecięcia stwarzając jednak możliwość odrzucenia potencjalnie dobrych punktów. Dobieramy jednak wartość $ eps $ jako liczbę bardzo małą, np. $ 0.00001 $ minimalizując szansę na odrzucenie poprawnego punktu przeciącia.

5
Twoja ocena: Brak Ocena: 5 (2 ocen)

Copyright © 2008-2010 Wrocławski Portal Informatyczny

design: rafalpolito.com