C++ DOT函数代码示例

/* compute oriented extents of a ellipsoid */
void ELLIP_Oriented_Extents (ELLIP *eli, double *vx, double *vy, double *vz, double *extents)
{
  double *c = eli->cur_center,
	 *r = eli->cur_sca,
	 *R = eli->cur_rot,
          p [8][3] = {{- r[0], - r[1], - r[2]},
                      {+ r[0], - r[1], - r[2]}, 
                      {+ r[0], + r[1], - r[2]}, 
                      {- r[0], + r[1], - r[2]}, 
                      {- r[0], - r[1], + r[2]},
                      {+ r[0], - r[1], + r[2]}, 
                      {+ r[0], + r[1], + r[2]}, 
                      {- r[0], + r[1], + r[2]}},
	  q [3],
	  e [3];

  for (int i = 0; i < 8; i ++)
  {
    NVADDMUL (c, R, p[i], q);

    e [0] = DOT (vx, q);
    e [1] = DOT (vy, q);
    e [2] = DOT (vz, q);

    if (e [0] < extents [0]) extents [0] = e [0];
    if (e [1] < extents [1]) extents [1] = e [1];
    if (e [2] < extents [2]) extents [2] = e [2];
    if (e [0] > extents [3]) extents [3] = e [0];
    if (e [1] > extents [4]) extents [4] = e [1];
    if (e [2] > extents [5]) extents [5] = e [2];
  }
}

static int intersect (unsigned int prim_id , const ray_t *ray, intersect_t *isect) {
    vector_t dest;
    double B, C, D;
    const float EPSILON = 0.003;

    SUB (dest, ray->orig, get_sphere_center (prim_id));
    B = DOT (dest, ray->dir);
    C = DOT (dest, dest);
    C -= get_sphere_radius (prim_id) * get_sphere_radius (prim_id);

    D = B*B  - C;

    if (D > 0) {
        const double sqrt_D = sqrt (D),
                     t0 = -B - sqrt_D,
                     t1 = -B + sqrt_D;

        if (t0 > EPSILON && t0 <= isect->t) {
            isect->t = t0;
            return 1;
        }
        else if (t1 > EPSILON && t1 <= isect->t) {
            isect->t = t1;
            return 1;
        }
        return 0;
    }
    return 0;
}

int
o_face(		/* determine if face intersects cube */
	OBJREC  *o,
	CUBE  *cu
)
{
	FVECT  cumin, cumax;
	FVECT  v1, v2;
	double  d1, d2;
	int  vloc;
	register FACE  *f;
	register int  i, j;
				/* get face arguments */
	f = getface(o);
	if (f->area == 0.0)		/* empty face */
		return(O_MISS);
					/* compute cube boundaries */
	for (j = 0; j < 3; j++)
		cumax[j] = (cumin[j] = cu->cuorg[j]-FTINY)
				+ cu->cusize + 2.0*FTINY;

	vloc = ABOVE | BELOW;		/* check vertices */
	for (i = 0; i < f->nv; i++)
		if ( (j = plocate(VERTEX(f,i), cumin, cumax)) )
			vloc &= j;
		else
			return(O_HIT);	/* vertex inside */

	if (vloc)			/* all to one side */
		return(O_MISS);
	
	for (i = 0; i < f->nv; i++) {	/* check edges */
		if ((j = i + 1) >= f->nv)
			j = 0;			/* wrap around */
		VCOPY(v1, VERTEX(f,i));		/* clip modifies */
		VCOPY(v2, VERTEX(f,j));		/* the vertices! */
		if (clip(v1, v2, cumin, cumax))
			return(O_HIT);		/* edge inside */
	}
					/* see if cube cuts plane */
	for (j = 0; j < 3; j++)
		if (f->norm[j] > 0.0) {
			v1[j] = cumin[j];
			v2[j] = cumax[j];
		} else {
			v1[j] = cumax[j];
			v2[j] = cumin[j];
		}
	if ((d1 = DOT(v1, f->norm) - f->offset) > FTINY)
		return(O_MISS);
	if ((d2 = DOT(v2, f->norm) - f->offset) < -FTINY)
		return(O_MISS);
					/* intersect face */
	for (j = 0; j < 3; j++)
		v1[j] = (v1[j]*d2 - v2[j]*d1)/(d2 - d1);
	if (inface(v1, f))
		return(O_HIT);
	
	return(O_MISS);		/* no intersection */
}

/* Update L-BFGS approximation of the Hessian. */
static void
update(opk_vmlmn_t* opt,
       const opk_vector_t* x,
       const opk_vector_t* g)
{
  double sty, yty;
  opk_index_t k;

  k = slot(opt, 0);
  AXPBY(S(k), 1, x, -1, opt->x0);
  AXPBY(Y(k), 1, g, -1, opt->g0);
  if (opt->method != OPK_VMLMN) {
    /* Compute initial inverse Hessian approximation. */
    sty = DOT(Y(k), S(k));
    if (sty <= 0) {
      RHO(k) = 0;
    } else {
      RHO(k) = 1/sty;
      yty = DOT(Y(k), Y(k));
      if (yty > 0) {
        opt->gamma = sty/yty;
      }
    }
  }
  ++opt->updates;
  if (opt->mp < opt->m) {
    ++opt->mp;
  }
}

void block::TBOMB2() {
	e[0] = e[4] = true;

	BOMB2(12);
	DOT(205, 8, 12, 12, 14);
	DOT(203, 8, 12, 13);
}

/* compute oriented extents of a sphere */
void SPHERE_Oriented_Extents (SPHERE *sph, double *vx, double *vy, double *vz, double *extents)
{
  double e [6], len [3], r;

  extents [0] = extents [1] = extents [2] =  DBL_MAX;
  extents [3] = extents [4] = extents [5] = -DBL_MAX;

  len [0] = LEN (vx);
  len [1] = LEN (vy);
  len [2] = LEN (vz);
    
  e [0] = DOT (sph->cur_center, vx);
  e [1] = DOT (sph->cur_center, vy);
  e [2] = DOT (sph->cur_center, vz);
  COPY (e, e + 3);
  r = sph->cur_radius;
  e [0] -= r / len [0];
  e [1] -= r / len [1];
  e [2] -= r / len [2];
  e [3] += r / len [0];
  e [4] += r / len [1];
  e [5] += r / len [2];

  if (e [0] < extents [0]) extents [0] = e [0];
  if (e [1] < extents [1]) extents [1] = e [1];
  if (e [2] < extents [2]) extents [2] = e [2];
  if (e [3] > extents [3]) extents [3] = e [3];
  if (e [4] > extents [4]) extents [4] = e [4];
  if (e [5] > extents [5]) extents [5] = e [5];
}

/* Update L-BFGS approximation of the Hessian. */
static void
update(opk_vmlmb_t* opt,
       const opk_vector_t* x,
       const opk_vector_t* g)
{
  double sty, yty;
  opk_index_t k;

  k = SLOT(0);
  AXPBY(S(k), 1, x, -1, opt->x0);
  AXPBY(Y(k), 1, g, -1, opt->g0);
  if (opt->method != OPK_VMLMB) {
    /* Compute initial inverse Hessian approximation. */
    sty = DOT(S(k), Y(k));
    if (sty <= 0) {
      /* This pair will be skipped.  This may however indicate a problem, see
         Nocedal & Wright "Numerical Optimization", section 8.1, p. 201 (1999).
         FIXME: restart? */
      RHO(k) = 0;
    } else {
      /* Compute RHO(k) and GAMMA. */
      RHO(k) = 1/sty;
      yty = DOT(Y(k), Y(k));
      if (yty > 0) {
        opt->gamma = sty/yty;
      }
    }
  }
  ++opt->updates;
  if (opt->mp < opt->m) {
    ++opt->mp;
  }
}

int Sphere::Intersect( Ray& a_Ray, float& a_Dist )
{
	vector3 v = a_Ray.GetOrigin() - m_Centre;
	float b = -DOT( v, a_Ray.GetDirection() );
	float det = (b * b) - DOT( v, v ) + m_SqRadius;
	int retval = MISS;
	if (det > 0)
	{
		det = sqrtf( det );
		float i1 = b - det;
		float i2 = b + det;
		if (i2 > 0)
		{
			if (i1 < 0) 
			{
				if (i2 < a_Dist) 
				{
					a_Dist = i2;
					retval = INPRIM;
				}
			}
			else
			{
				if (i1 < a_Dist)
				{
					a_Dist = i1;
					retval = HIT;
				}
			}
		}
	}
	return retval;
}

static inline bool8 AtiPlaneBoxOverlap( const TBM_FLOAT normal[3], const TBM_FLOAT d, const TBM_FLOAT maxbox[3])
{
   TBM_FLOAT vmin[3];
   TBM_FLOAT vmax[3];
   for (int32 q = X; q <= Z; q++)
   {
      if (normal[q] > TBM_ZERO)
      {
         vmin[q] = -maxbox[q];
         vmax[q] = maxbox[q];
      }
      else
      {
         vmin[q] = maxbox[q];
         vmax[q] = -maxbox[q];
      }
   }
   if (DOT (normal, vmin) + d > TBM_ZERO)
   {
      return FALSE;
   }
   if (DOT (normal, vmax) + d >= TBM_ZERO)
   {
      return TRUE;
   }
   return FALSE;
}

/* push 'p' deeper inside of convices bounded by two plane sets */
static int refine_point (double *pa, int npa, double *pb, int npb, double *p, double *epsout)
{
  double *pla, *end, eps, d, q [3];
  short pushed, iter, imax;

  imax = 4;
  iter = 0;
  eps = GEOMETRIC_EPSILON * (1 << (imax + 1));
  do
  {
    for (pushed = 0, pla = pa, end = pa + npa * 6; pla < end; pla += 6)
    {
      SUB (p, pla + 3, q);
      d = DOT (pla, q);
      if (d > -GEOMETRIC_EPSILON) { SUBMUL (p, eps, pla, p); pushed = 1; }
    }

    for (pla = pb, end = pb + npb * 6; pla < end; pla += 6)
    {
      SUB (p, pla + 3, q);
      d = DOT (pla, q);
      if (d > -GEOMETRIC_EPSILON) { SUBMUL (p, eps, pla, p); pushed = 1; }
    }

    eps *= 0.5;

  } while (pushed && iter ++ < imax);

  *epsout = 10 * eps;

  return !pushed;
}

/* Find convex hull vertex to complete triangle (oriented call) */
static RBFNODE *
find_chull_vert(const RBFNODE *rbf0, const RBFNODE *rbf1)
{
	FVECT	vmid, vejn, vp;
	RBFNODE	*rbf, *rbfbest = NULL;
	double	dprod, area2, bestarea2 = FHUGE, bestdprod = -.5;

	VSUB(vejn, rbf1->invec, rbf0->invec);
	VADD(vmid, rbf0->invec, rbf1->invec);
	if (normalize(vejn) == 0 || normalize(vmid) == 0)
		return(NULL);
						/* XXX exhaustive search */
	/* Find triangle with minimum rotation from perpendicular */
	for (rbf = dsf_list; rbf != NULL; rbf = rbf->next) {
		if ((rbf == rbf0) | (rbf == rbf1))
			continue;
		tri_orient(vp, rbf0->invec, rbf1->invec, rbf->invec);
		if (DOT(vp, vmid) <= FTINY)
			continue;		/* wrong orientation */
		area2 = .25*DOT(vp,vp);
		VSUB(vp, rbf->invec, vmid);
		dprod = -DOT(vp, vejn);
		VSUM(vp, vp, vejn, dprod);	/* above guarantees non-zero */
		dprod = DOT(vp, vmid) / VLEN(vp);
		if (dprod <= bestdprod + FTINY*(1 - 2*(area2 < bestarea2)))
			continue;		/* found better already */
		if (overlaps_tri(rbf0, rbf1, rbf))
			continue;		/* overlaps another triangle */
		rbfbest = rbf;
		bestdprod = dprod;		/* new one to beat */
		bestarea2 = area2;
	}
	return(rbfbest);
}

// Super fire item (fogo level 9)
void block::SFIREIT() {
	e[0] = e[3] = true;
	item = 'F';

	CIRCLE(12, 14);
	DOT(NR, 14, 0, 11, 15);
	DOT(DR, 12, 14, 12, 14);
}

// Time bomb item
void block::TBOMBIT() {
	e[0] = e[3] = true;
	item = 't';

	BOMB1(14);
	DOT(205, 8, 1, 12, 14);
	DOT(203, 8, 1, 13);
}

// Forma de cruz
void block::FIRECENTER() {
	e[0] = e[7] = true;

	BLOCK(NR, 15, 0);
	DOT(DR, 14, 12, 11, 15);
	DOT(NR, 14, 0, 12, 14, 32, 34);
	DOT(UR, 14, 12, 31, 35);
}

void Transform(MATRIX a, VECTOR u, VECTOR v)
{
  /* Applies matrix "a" to vector "u", returning vector "v" */

  v[X] = DOT(a[X], u);
  v[Y] = DOT(a[Y], u);
  v[Z] = DOT(a[Z], u);
}