მელერ-ტრამბორის კვეთის ალგორითმი

**მელერ-ტრამბორის ალგორითმით** დგინდება სხივის კვეთა სამკუთხედთან სამ განზომილებაში. თომას მელერისა და ბენ ტრამბორის ეს ალგორითმი არ საჭიროებს წერტილების სიბრტყეზე გადატანას, ზოგავს მეხსიერებას და თავის თანამედროვე იმპლემენტაციებთან შედარებით სწრაფია. ის პუბლიკაცია ACM-მა გამოსცა 1997 წელს სახელით „სწრაფი, მინიმალური მეხსიერების სამკუთხედ-სხივის კვეთა“ და კომპიუტერული გრაფიკის სფეროში ერთ-ერთ ყველაზე ცნობილ ნაშრომს წარმოადგენს. იმის მიუხედავად, რომ უკეთესი ალგორითმებიც არსებობს, მელერ-ტრამბორის ალგორითმი დღესაც ფართოდაა გავრცელებული მისი სიმარტივის გამო.

## ალგორითმი 
### წერტილი სხივზე: 
R(t) = O + tD (1)

O - საწყისი წერილია, t - მანძილი, D - სხივის მიმართულება

### წერტილი სამკუთხედზე: 
T(u,v) = (1 - u -v)V_0 + uV_1 + vV_2 (2)

V_0,V_1,V_2 სამკუთხედის წვეროს კოორდინატებია

### მთავარი ტოლობა: 
სხივი რომ სამკუთხედს კვეთდეს უნდა არსებობდეს წერტილი სხივზე, რომელიც სამკუთხედს ეკუთვნის.

R(t) = T(u,v) ანუ:

O + tD = (1 - u - v)V_0 + uV_1 + vV_2 (3)

გადავანაცვლოთ ცვლადები და ვიღებთ:

[ -D, V_1 - V_0, V_2 - V_0]\begin{bmatrix} t \\ u \\ v \end{bmatrix} = O - V_0 (4)

კვეთის სანახავად უნდა ამოიხსნას განტოლება (4).

### ამოხსნა: 
E_1 = V_1 - V_0

E_2 = V_2 - V_0

T = O - V_0

P = (D\times E_2)

Q = T \times E_1

განტოლების (4) ამოხსნა ხდება შემდეგნაირად:

\begin{bmatrix} t \\ u \\ v \end{bmatrix} =
\frac{1}{(D \times E_2) \cdot E_1} 
\begin{bmatrix} 
(T \times E_1)\cdot E_2 \\
(D \times E_2)\cdot T  \\
(T \times E_1)\cdot D 
\end{bmatrix} =
\frac{1}{P \times E_1}
\begin{bmatrix} 
(Q \times E_2)\\
(P \times T)\\
(Q \times D)
\end{bmatrix}

 (6)

მივიღებთ სამეულს (t, u, v)

## იმპლემენტაცია 
```c
/* Ray-Triangle Intersection Test Routines          */
/* Different optimizations of my and Ben Trumbore's */
/* code from journals of graphics tools (JGT)       */
/* http://www.acm.org/jgt/                          */
/* by Tomas Moller, May 2000                        */

#include 

#define EPSILON 0.000001
#define CROSS(dest,v1,v2) \
          dest[0]=v1[1]*v2[2]-v1[2]*v2[1]; \
          dest[1]=v1[2]*v2[0]-v1[0]*v2[2]; \
          dest[2]=v1[0]*v2[1]-v1[1]*v2[0];
#define DOT(v1,v2) (v1[0]*v2[0]+v1[1]*v2[1]+v1[2]*v2[2])
#define SUB(dest,v1,v2) \
          dest[0]=v1[0]-v2[0]; \
          dest[1]=v1[1]-v2[1]; \
          dest[2]=v1[2]-v2[2]; 

/* the original jgt code */
int intersect_triangle(double orig[3], double dir[3],
		       double vert0[3], double vert1[3], double vert2[3],
		       double *t, double *u, double *v)
{
   double edge1[3], edge2[3], tvec[3], pvec[3], qvec[3];
   double det,inv_det;

   /* find vectors for two edges sharing vert0 */
   SUB(edge1, vert1, vert0);
   SUB(edge2, vert2, vert0);

   /* begin calculating determinant - also used to calculate U parameter */
   CROSS(pvec, dir, edge2);

   /* if determinant is near zero, ray lies in plane of triangle */
   det = DOT(edge1, pvec);

   if (det > -EPSILON && det  1.0)
     return 0;

   /* prepare to test V parameter */
   CROSS(qvec, tvec, edge1);

   /* calculate V parameter and test bounds */
   *v = DOT(dir, qvec) * inv_det;
   if (*v  1.0)
     return 0;

   /* calculate t, ray intersects triangle */
   *t = DOT(edge2, qvec) * inv_det;

   return 1;
}

/* code rewritten to do tests on the sign of the determinant */
/* the division is at the end in the code                    */
int intersect_triangle1(double orig[3], double dir[3],
			double vert0[3], double vert1[3], double vert2[3],
			double *t, double *u, double *v)
{
   double edge1[3], edge2[3], tvec[3], pvec[3], qvec[3];
   double det,inv_det;

   /* find vectors for two edges sharing vert0 */
   SUB(edge1, vert1, vert0);
   SUB(edge2, vert2, vert0);

   /* begin calculating determinant - also used to calculate U parameter */
   CROSS(pvec, dir, edge2);

   /* if determinant is near zero, ray lies in plane of triangle */
   det = DOT(edge1, pvec);

   if (det > EPSILON)
   {
      /* calculate distance from vert0 to ray origin */
      SUB(tvec, orig, vert0);
      
      /* calculate U parameter and test bounds */
      *u = DOT(tvec, pvec);
      if (*u  det)
	 return 0;
      
      /* prepare to test V parameter */
      CROSS(qvec, tvec, edge1);
      
      /* calculate V parameter and test bounds */
      *v = DOT(dir, qvec);
      if (*v  det)
	 return 0;
      
   }
   else if(det  0.0 || *u  0.0 || *u + *v  EPSILON)
   {
      /* calculate U parameter and test bounds */
      *u = DOT(tvec, pvec);
      if (*u  det)
	 return 0;
      
      /* prepare to test V parameter */
      CROSS(qvec, tvec, edge1);
      
      /* calculate V parameter and test bounds */
      *v = DOT(dir, qvec);
      if (*v  det)
	 return 0;
      
   }
   else if(det  0.0 || *u  0.0 || *u + *v  EPSILON)
   {
      *u = DOT(tvec, pvec);
      if (*u  det)
	 return 0;
            
      /* calculate V parameter and test bounds */
      *v = DOT(dir, qvec);
      if (*v  det)
	 return 0;
      
   }
   else if(det  0.0 || *u  0.0 || *u + *v < det)
	 return 0;
   }
   else return 0;  /* ray is parallell to the plane of the triangle */

   *t = DOT(edge2, qvec) * inv_det;
   (*u) *= inv_det;
   (*v) *= inv_det;

   return 1;
}
```

## სქოლიო 
კატეგორია:გეომეტრიული ალგორითმები