C++ Vec3f::Dot3方法代码示例

bool Face::plane_intersect(const Ray &r, Hit &h, bool intersect_backfacing, bool* backfacing_hit) {

	// insert the explicit equation for the ray into the implicit equation of the plane

	// equation for a plane
	// ax + by + cz = d;
	// normal . p + direction = 0
	// plug in ray
	// origin + direction * t = p(t)
	// origin . normal + t * direction . normal = d;
	// t = d - origin.normal / direction.normal;

	Vec3f normal = computeNormal();
	double d = normal.Dot3((*this)[0]->get());

	double numer = d - r.getOrigin().Dot3(normal);
	double denom = r.getDirection().Dot3(normal);

	if (denom == 0) return 0;	// parallel to plane

	if (!intersect_backfacing && normal.Dot3(r.getDirection()) >= 0) 
		return 0; // hit the backside

	double t = numer / denom;
	if (t > EPSILON && t < h.getT()) {
		h.set(t,this->getMaterial(),normal,this);
		assert (h.getT() >= EPSILON);
		//hit the backside but that's okay in this case
		if (normal.Dot3(r.getDirection()) >= 0){
			*backfacing_hit = true;
		}
		return 1;
	}
	return 0;
}

bool Plane::intersect(const Ray & r, Hit & h, float tmin)
{
	assert (tmin >= 0.0f);

	Vec3f org = r.getOrigin();
	Vec3f dir = r.getDirection();

	if (dir.Dot3(normal) < 1.0e-7 && dir.Dot3(normal) > -1.0e-7) {
		return false;
	}					// Appromately parrell to plane

	float tempT = (offset - org.Dot3(normal)) / dir.Dot3(normal);
	if (tempT <= 1e-6) {
		return false;
	}
	else if (tempT >= tmin && tempT < h.getT()) 
	{
		// Update Hit Point
		normal.Normalize();
		h.set(tempT, NULL, normal, color, r);

		return true;
	}
	
	return false;
}

bool Sphere::intersect(const Ray &r, Hit &h, float tmin)
{
	Vec3f v = center - r.getOrigin();
	float tp = v.Dot3(r.getDirection());
	float det = tp*tp - v.Dot3(v) + radius*radius;
	//intersect
	if(det > 0)
	{
		//t'
		det = sqrtf(det); 

		float t1 = tp - det;
		float t2 = tp + det;

		if(t1 > tmin && t1 < h.getT())
		{
			Vec3f normal = (r.pointAtParameter(t1) - center);
			normal /= radius;
			normal.Normalize();
			h.set(t1,material,normal,r);
			return 1;
		}
		else if(t2 > tmin && t2 < h.getT())
		{
			//sphere's normal
			Vec3f normal = (r.pointAtParameter(t2) - center);
			normal /= radius;
			normal.Normalize();
			h.set(t2,material,normal,r);
			return 1;
		}
	}
	return 0;
}

float sphere::testIntersection(Vec3f eye, Vec3f dir)
{
	float disc = pow(dir.Dot3(eye - center), 2) - dir.Dot3(dir) * ((eye - center).Dot3(eye - center) - radius * radius);
	if(disc < 0) return MAX_DIST;
	float t1 = (-1 * dir.Dot3(eye - center) + sqrt(disc)) / dir.Dot3(dir);
	float t2 = (-1 * dir.Dot3(eye - center) - sqrt(disc)) / dir.Dot3(dir);
	if(t1 < 0 || t2 < 0)
		return MAX_DIST;
	return min(t1, t2);
}

float Edge::DihedralAngle(){
  /*
   * Warning this function returns NULL 
   * when there is an edge without two 
   * adjcent triangles
   *
   */

  // Is there even an angle here?
  if(opposite == NULL)
    return (float)NULL;

  // Find the angle of the face of a triangle
  Triangle* a = triangle;
  Triangle* b = opposite->getTriangle();

  Vec3f normalA = a->getNormal();
  Vec3f normalB = b->getNormal();

  // Using Equation theta = acos( (a . b) / (|a||b|)) 
  double top = normalA.Dot3(normalB);
  double bottom = normalA.Length() * normalB.Length();
  double result = acos(top/bottom);
  return result;

}

Vec3f mirrorDirection(const Vec3f &normal, const Vec3f &incoming)
{
	//(反射光方向：R = V - 2(V.N)N )
	Vec3f reflectDir = incoming - 2 * (incoming.Dot3(normal))*normal;
	reflectDir.Normalize();
	return reflectDir;
}

Vec3f RayTracer::shadow(const Vec3f &point,
			const Vec3f &pointOnLight,
			const Face *f,
			const Ray &ray,
			const Hit &hit) const
{
        const Vec3f normal(hit.getNormal());
        const Material *m = hit.getMaterial();

	Vec3f dirToLight = pointOnLight - point;
	dirToLight.Normalize();
	/* If dot product < 0, surface is not facing light */
	if (normal.Dot3(dirToLight) > 0) {
		Ray rayToLight(point, dirToLight);
		Hit hLight;
		bool blocked = CastRay(rayToLight, hLight, false);
		while (std::fabs(hLight.getT()) < SURFACE_EPSILON &&
				std::fabs((pointOnLight - point).Length()) > SURFACE_EPSILON) {
			rayToLight = Ray(rayToLight.pointAtParameter(SURFACE_EPSILON),
					dirToLight);
			blocked = CastRay(rayToLight, hLight, false);
		}
		if (hLight.getT() == FLT_MAX || hLight.getMaterial() != f->getMaterial()) {
			return Vec3f(0, 0, 0);
		}

		const Vec3f lightColor = 0.2 * f->getMaterial()->getEmittedColor() * f->getArea();
		return m->Shade(ray,hit,dirToLight,lightColor,args);
	}
	return Vec3f(0, 0, 0);
}

Vec3f RayTracer::mirrorDirection(const Vec3f &normal, const Vec3f &incoming) const
{
	Vec3f reflectionDir;
	reflectionDir = incoming - normal * (incoming.Dot3(normal)) * 2;
	reflectionDir.Normalize();
	return reflectionDir;
}

bool Sphere::intersect(const Ray &r, Hit &h, float tmin)
{
	//直线方程 P = t*D + R    ①
	//圆方程  ||P||= raduis   ②
	//将①带入②, 由点乘的性质(满足分配率，交换律)，求得t^2+2RDt+R^2-r^2=0
	//得t = -RD±sqrt(RD^2-R^2+r^2)
	//选择距离较近的那个点t = -RD-sqrt(RD^2-R^2+r^2)
	Vec3f D = r.getDirection();
	Vec3f R = r.getOrigin()-center;


	float R2 = R.Dot3(R);
	float RD = R.Dot3(D);
	
	float b2_4ac = RD*RD-R2+radius*radius;//R2 RD RDRD这些其实可以存在ray里，但是算法研究以程序清晰为要
	if(b2_4ac<0)return 0;
	float t;
	t = -RD - sqrt(b2_4ac);
	if(t<0){
		t = -RD + sqrt(b2_4ac);
		if(t < 0)return 0;
	}
	if(t<h.getT())
	{
		h.set(t, &color, r);
	}

	return 1;
}

void PhotonMapping::TracePhoton(const Vec3f &position, const Vec3f &direction, 
				const Vec3f &energy, int iter) {
  
  if(iter>args->num_bounces){
    return;
  }

  Hit h = Hit();
  Ray R = Ray(position, direction);
  bool intersect = raytracer->CastRay(R, h, false);
  if(!intersect){
    return;
  }
  Material *m = h.getMaterial();
  Vec3f normal = h.getNormal();
  Vec3f point = R.pointAtParameter(h.getT());
  Vec3f opDirec = direction;
  opDirec.Negate();
  opDirec.Normalize();
  Vec3f diffuse = m->getDiffuseColor(), reflec = m->getReflectiveColor();
  double diffuseAnswer = diffuse.x()+diffuse.y()+diffuse.z();
  double reflecAnswer = reflec.x()+reflec.y()+reflec.z();
  double total = reflecAnswer+diffuseAnswer;
  diffuseAnswer /= total;
  reflecAnswer /= total;
  double seed = GLOBAL_mtrand.rand();
  if(seed <= diffuseAnswer && seed >= 0){
    Vec3f newEnergy = energy * diffuse;
    Vec3f newPosition = point;
    Vec3f newDirection = Vec3f(GLOBAL_mtrand.rand(),GLOBAL_mtrand.rand(),GLOBAL_mtrand.rand());
    newDirection.Normalize();
    Photon answer = Photon(point,opDirec,newEnergy,iter+1);
    kdtree->AddPhoton(answer);
    TracePhoton(newPosition, newDirection, newEnergy, iter+1);
  }
  else if(seed>diffuseAnswer && seed <= 1){
    Vec3f newEnergy = energy * reflec;
    Vec3f newPosition = point;
    Vec3f newDirection = direction - 2 * direction.Dot3(normal) * normal;
    Photon answer = Photon(point,opDirec,newEnergy,iter+1);
    kdtree->AddPhoton(answer);
    TracePhoton(newPosition, newDirection, newEnergy, iter+1);
  }
  // ==============================================
  // ASSIGNMENT: IMPLEMENT RECURSIVE PHOTON TRACING
  // ==============================================

  // Trace the photon through the scene.  At each diffuse or
  // reflective bounce, store the photon in the kd tree.

  // One optimization is to *not* store the first bounce, since that
  // direct light can be efficiently computed using classic ray
  // tracing.



}

Vec3f PhotonMapping::GatherIndirect(const Vec3f &point, const Vec3f &normal, const Vec3f &direction_from) const {


  if (kdtree == NULL) { 
    std::cout << "WARNING: Photons have not been traced throughout the scene." << std::endl;
    return Vec3f(0,0,0); 
  }


  // ================================================================
  // ASSIGNMENT: GATHER THE INDIRECT ILLUMINATION FROM THE PHOTON MAP
  // ================================================================


  // collect the closest args->num_photons_to_collect photons
  // determine the radius that was necessary to collect that many photons
  // average the energy of those photons over that radius
  double radius=1;

  std::vector<Photon> potentials;
  while (potentials.size()<args->num_photons_to_collect)
  {
	  BoundingBox b(point-Vec3f(radius,radius,radius),point+Vec3f(radius,radius,radius));

	  kdtree->CollectPhotonsInBox(b,potentials);
	  for (int i=0;i<potentials.size();i++)
	  {
		  if (normal.Dot3(potentials[i].getDirectionFrom())>0)
		  {
			  potentials.erase(potentials.begin()+i);
			  i--;
		  }
		  else if ((potentials[i].getPosition()-point).Length()>radius)
		  {
			  potentials.erase(potentials.begin()+i);
			  i--;
		  }
	  }

	  radius*=2;


  }
  assert(potentials.size()>=args->num_photons_to_collect);
  for (unsigned int i=0;i<potentials.size();i++)
	  potentials[i].setPoint(point);
  std::sort(potentials.begin(),potentials.end(),comparePhotons);

  radius = (potentials[args->num_photons_to_collect-1].getPosition()-point).Length()/2;
  Vec3f color;
  for (int i=0;i<args->num_photons_to_collect;i++)
	  color+=potentials[i].getEnergy();
  color*=1/(radius*radius);
  //std::cout<<radius<<"\n";
  // return the color
  return color;
}

Vec3f PhongMaterial::Shade(const Ray& ray, const Hit& hit, const Vec3f& dirToLight, const Vec3f& lightColor) const {
	Vec3f v = ray.getDirection() * -1.0f;
    Vec3f n = hit.getNormal();
	// Vec3f h = dirToLight + ray.getDirection() * -1.0f;
	// h.Normalize();
	float diffuse = dirToLight.Dot3(n);
	Vec3f r = n * 2.0 * n.Dot3(dirToLight) - dirToLight;
	float specular = v.Dot3(r);
	
	if (diffuse < 0.0f) diffuse = 0.0f;
	if (specular < 0.0f) specular = 0.0f;
	
	Vec3f color = this->getDiffuseColor();
	color = color * diffuse + m_specularColor * pow(specular, m_exponent);
	color.Scale(lightColor);
	
	return color;
}

float triangle::testIntersection(Vec3f eye, Vec3f dir)
{
	//see the book/slides for a description of how to use Cramer's rule to solve
	//for the intersection(s) of a line and a plane, implement it here and
	//return the minimum distance (if barycentric coordinates indicate it hit
	//the triangle) otherwise 9999999
	
	Vec3f edge1 =  alpha;
	Vec3f edge2 =  beta;

	Vec3f pvec;
	Vec3f::Cross3(pvec,dir, edge2);
	
	float det = edge1.Dot3(pvec);
	
	// no intersection
	if(det > (-1 * EPSILON) && det < EPSILON)
            return MISS;
	
	float invDet = 1/det;
	Vec3f tvec = eye-point0;
	
	float u = tvec.Dot3(pvec) * invDet;
	
	if(u < 0 || u > 1)
		return MISS;
	
	Vec3f qvec;
	Vec3f::Cross3(qvec,tvec, edge1);
        
	float v = dir.Dot3(qvec) * invDet;
	
	if(v < 0 || (u + v) > 1)
		return MISS;          
      
	float dist = edge2.Dot3(qvec) * invDet;
	
	if(dist < 0)
		return MISS;
	
	return dist;
}

bool Sphere::intersect(const Ray &r, Hit &h, float tmin)
{
    Vec3f temp = r.getOrigin() - mCenterPoint;
    Vec3f rayDirection = r.getDirection();

    double a = rayDirection.Dot3(rayDirection);
    double b = 2*rayDirection.Dot3(temp);
    double c = temp.Dot3(temp) - mRadius*mRadius;

#ifdef DEBUG
    //cout << "temp="  << temp<<endl;
    //printf("Sphere::intersect, a=%f, b=%f, c=%f\n", a, b, c);
#endif

    double discriminant = b*b - 4*a*c;

    if (discriminant > 0)
    {
        discriminant = sqrt(discriminant);
        double t = (- b - discriminant) / (2*a);

        if (t < tmin)
            t = (- b + discriminant) / (2*a);

        if (t < tmin || t > T_MAX)
            return false;

#ifdef DEBUG
      //  printf("Sphere::intersect, there is a hit, t=%f\n", t);
#endif
        Vec3f normal = r.getOrigin() + t * r.getDirection() - mCenterPoint;
        normal.Normalize();
        h.set(t, mMaterial, normal, r);
        return true;
    }

    return false;
}

//何时返回true？非全反射时 何时返回false
bool RayTracer::transmittedDirection(const Vec3f &normal, const Vec3f &incoming, float index_i, float index_t, Vec3f &transmitted) const
{
	float nr = index_i / index_t;
	Vec3f I = incoming*(-1.0f);
	float cosI = I.Dot3(normal);
	float isAllTrans = 1 - nr*nr*(1 - cosI*cosI);
	if(isAllTrans < 0)  //全反射
		return false;

	float cosT = sqrt(isAllTrans);
	transmitted = normal*(nr*cosI - cosT) - I*nr;
	transmitted.Normalize();
	return true;
}