/* Definitions et fonctions d'utilite generale  */
/* pour la gestion de directions                */
/* dans un espace a 3 dimensions                */
/*                                              */
/* Auteur: Nicolas JANEY                        */
/* nicolas.janey@univ-fcomte.fr                 */
/* Novembre 2009                                */

#include <math.h>

#include "RayTracing-Direction.h"

/* Fonction de calcul de la direction du rayon  */
/* lumineux reflechi cree par un rayon lumineux */
/* incident defini par sa direction             */
/* de propagation                               */

void reflexion(direction *r,direction *i,direction *n) {
  double ps2 = 2.0*produitScalaire(n,i);
  r->x = ps2*n->x - i->x;
  r->y = ps2*n->y - i->y;
  r->z = ps2*n->z - i->z;
  r->t = 0.0;
}

/* Fonction de calcul de la direction du rayon  */
/* lumineux transmis eventuellement cree        */
/* par un rayon lumineux incident defini        */
/* par sa direction de propagation              */

int transmission(direction *t,direction *i,direction *n,
                 double niSurnt) {
  double ps = produitScalaire(n,i);
  double v = 1.0-(niSurnt*niSurnt*(1.0-ps*ps));
  if ( v < 0.0 )
    return(0);
  v = niSurnt*ps-sqrt(v);
  t->x = v*n->x - niSurnt*i->x;
  t->y = v*n->y - niSurnt*i->y;
  t->z = v*n->z - niSurnt*i->z;
  t->t = 0.0;
  return(1);
}

/* Fonction de normalisation d'une direction    */

void normalise(direction *d) {
  double l = sqrt(d->x*d->x+d->y*d->y+d->z*d->z);
  if ( l != 0 ) {
    d->x /= l;
    d->y /= l;
    d->z /= l; }
}

/* Fonction de calcul du produit scalaire       */
/* entre deux directions                        */

double produitScalaire(direction *d1,direction *d2) {
  return(d1->x*d2->x+d1->y*d2->y+d1->z*d2->z);
}

/* Fonction de calcul de la direction           */
/* normalisee definie d'une position            */
/* a une autre position                         */

void calculDirectionNormalisee(direction *d,position *p1,position *p2) {
  d->x = p2->x-p1->x;
  d->y = p2->y-p1->y;
  d->z = p2->z-p1->z;
  d->t = 0.0;
  normalise(d);
}