C++ Intersection::GetBSDF方法代码示例

Spectrum DiffusePRTIntegrator::Li(const Scene *scene, const Renderer *,
            const RayDifferential &ray, const Intersection &isect,
            const Sample *sample, MemoryArena &arena) const {
    Spectrum L = 0.f;
    Vector wo = -ray.d;
    // Compute emitted light if ray hit an area light source
    L += isect.Le(wo);

    // Evaluate BSDF at hit point
    BSDF *bsdf = isect.GetBSDF(ray, arena);
    const Point &p = bsdf->dgShading.p;
    const Normal &n = bsdf->dgShading.nn;
    // Compute reflected radiance using diffuse PRT

    // Project diffuse transfer function at point to SH
    Spectrum *c_transfer = arena.Alloc<Spectrum>(SHTerms(lmax));
    SHComputeDiffuseTransfer(p, Faceforward(n, wo), isect.rayEpsilon,
        scene, *sample->rng, nSamples, lmax, c_transfer);

    // Compute integral of product of incident radiance and transfer function
    Spectrum LT = 0.f;
    for (int i = 0; i < SHTerms(lmax); ++i)
        LT += c_in[i] * c_transfer[i];

    // Compute reflectance at point for diffuse transfer
    const int sqrtRhoSamples = 6;
    float rhoRSamples[2*sqrtRhoSamples*sqrtRhoSamples];
    StratifiedSample2D(rhoRSamples, sqrtRhoSamples, sqrtRhoSamples, *sample->rng);
    Spectrum Kd = bsdf->rho(wo, sqrtRhoSamples*sqrtRhoSamples, rhoRSamples,
        BSDF_ALL_REFLECTION) * INV_PI;
    return L + Kd * LT.Clamp();
}

static Spectrum L(const Scene *scene, const Renderer *renderer,
        const Camera *camera, MemoryArena &arena, RNG &rng, int maxDepth,
        bool ignoreDirect, const MLTSample &sample) {
    // Generate camera ray from Metropolis sample
    RayDifferential ray;
    float cameraWeight = camera->GenerateRayDifferential(sample.cameraSample,
                                                         &ray);
    Spectrum pathThroughput = cameraWeight, L = 0.;
    bool specularBounce = false, allSpecular = true;
    for (int pathLength = 0; pathLength < maxDepth; ++pathLength) {
        // Find next intersection in Metropolis light path
        Intersection isect;
        if (!scene->Intersect(ray, &isect)) {
            bool includeLe = ignoreDirect ? (specularBounce && !allSpecular) :
                                            (pathLength == 0 || specularBounce);
            if (includeLe)
                for (uint32_t i = 0; i < scene->lights.size(); ++i)
                   L += pathThroughput * scene->lights[i]->Le(ray);
            break;
        }
        if (ignoreDirect ? (specularBounce && !allSpecular) :
                           (specularBounce || pathLength == 0))
            L += pathThroughput * isect.Le(-ray.d);
        BSDF *bsdf = isect.GetBSDF(ray, arena);
        const Point &p = bsdf->dgShading.p;
        const Normal &n = bsdf->dgShading.nn;
        Vector wo = -ray.d;
        const PathSample &ps = sample.pathSamples[pathLength];
        // Sample direct illumination for Metropolis path vertex
        if (!ignoreDirect || pathLength > 0) {
            LightSample lightSample(ps.lightDir0, ps.lightDir1, ps.lightNum0);
            BSDFSample bsdfSample(ps.bsdfLightDir0, ps.bsdfLightDir1,
                                  ps.bsdfLightComponent);
            uint32_t lightNum = Floor2Int(ps.lightNum1 * scene->lights.size());
            lightNum = min(lightNum, (uint32_t)(scene->lights.size()-1));
            const Light *light = scene->lights[lightNum];
            L += pathThroughput *
                 EstimateDirect(scene, renderer, arena, light, p, n, wo,
                     isect.rayEpsilon, sample.cameraSample.time, bsdf, rng,
                     lightSample, bsdfSample);
        }

        // Sample direction for outgoing Metropolis path direction
        BSDFSample outgoingBSDFSample(ps.bsdfDir0, ps.bsdfDir1,
                                      ps.bsdfComponent);
        Vector wi;
        float pdf;
        BxDFType flags;
        Spectrum f = bsdf->Sample_f(wo, &wi, outgoingBSDFSample,
                                    &pdf, BSDF_ALL, &flags);
        if (f.IsBlack() || pdf == 0.)
            break;
        specularBounce = (flags & BSDF_SPECULAR) != 0;
        allSpecular &= specularBounce;
        pathThroughput *= f * AbsDot(wi, n) / pdf;
        ray = RayDifferential(p, wi, ray, isect.rayEpsilon);
        //pathThroughput *= renderer->Transmittance(scene, ray, NULL, rng, arena);
    }
    return L;
}

Spectrum DipoleSubsurfaceIntegrator::Li(const Scene *scene, const Renderer *renderer,
        const RayDifferential &ray, const Intersection &isect,
        const Sample *sample, RNG &rng, MemoryArena &arena) const {
    Spectrum L(0.);
    Vector wo = -ray.d;
    // Compute emitted light if ray hit an area light source
    L += isect.Le(wo);

    // Evaluate BSDF at hit point
    BSDF *bsdf = isect.GetBSDF(ray, arena);
    const Point &p = bsdf->dgShading.p;
    const Normal &n = bsdf->dgShading.nn;

    Spectrum rho_dr = Spectrum(1.0f);

    // Evaluate BSSRDF and possibly compute subsurface scattering
    BSSRDF *bssrdf = isect.GetBSSRDF(ray, arena);
    if (bssrdf && octree) {	
        Spectrum sigma_a  = bssrdf->sigma_a();
        Spectrum sigmap_s = bssrdf->sigma_prime_s();
        Spectrum sigmap_t = sigmap_s + sigma_a;
        if (!sigmap_t.IsBlack()) {
            // Use hierarchical integration to evaluate reflection from dipole model
            PBRT_SUBSURFACE_STARTED_OCTREE_LOOKUP(const_cast<Point *>(&p));
            DiffusionReflectance Rd(sigma_a, sigmap_s, bssrdf->eta());
            Spectrum Mo = octree->Mo(octreeBounds, p, Rd, maxError);
            FresnelDielectric fresnel(1.f, bssrdf->eta());
            Spectrum Ft = Spectrum(1.f) - fresnel.Evaluate(AbsDot(wo, n));
            float Fdt = 1.f - Fdr(bssrdf->eta());

	    // modulate SSS contribution by rho_dr
            //L += (INV_PI * Ft) * (Fdt * Mo);
	    rho_dr = wet->integrate_BRDF(bsdf, ray.d, 10, BxDFType(BSDF_REFLECTION | BSDF_GLOSSY));
	    L += (INV_PI * Ft) * (Fdt * Mo) * (Spectrum(1.0f) - rho_dr);
	    //L += (INV_PI * Ft) * (Fdt * Mo) * (Spectrum(0.0f));
	    
            PBRT_SUBSURFACE_FINISHED_OCTREE_LOOKUP();
        }
    }

    L += UniformSampleAllLights(scene, renderer, arena, p, n,
        wo, isect.rayEpsilon, ray.time, bsdf, sample, rng, lightSampleOffsets,
        bsdfSampleOffsets);

    if (ray.depth < maxSpecularDepth) {
        // Trace rays for specular reflection and refraction.

      //TODO: this has no effect?
        L += SpecularReflect(ray, bsdf, rng, isect, renderer, scene,
			     sample, arena);
        L += SpecularTransmit(ray, bsdf, rng, isect,
			      renderer, scene, sample, arena);
    }
    return L;
}

Spectrum IrradianceCacheIntegrator::pathL(Ray &r, const Scene *scene,
        const Renderer *renderer, const Sample *sample, MemoryArena &arena) const {
    Spectrum L(0.f);
    Spectrum pathThroughput = 1.;
    RayDifferential ray(r);
    bool specularBounce = false;
    for (int pathLength = 0; ; ++pathLength) {
        // Find next vertex of path
        Intersection isect;
        if (!scene->Intersect(ray, &isect))
            break;
        if (pathLength == 0)
            r.maxt = ray.maxt;
        else if (pathLength == 1)
            pathThroughput *= renderer->Transmittance(scene, ray, sample, arena, NULL);
        else
            pathThroughput *= renderer->Transmittance(scene, ray, NULL, arena, sample->rng);
        // Possibly add emitted light at path vertex
        if (specularBounce)
            L += pathThroughput * isect.Le(-ray.d);
        // Evaluate BSDF at hit point
        BSDF *bsdf = isect.GetBSDF(ray, arena);
        // Sample illumination from lights to find path contribution
        const Point &p = bsdf->dgShading.p;
        const Normal &n = bsdf->dgShading.nn;
        Vector wo = -ray.d;
        L += pathThroughput *
            UniformSampleOneLight(scene, renderer, arena, p, n, wo, isect.rayEpsilon,
                                  bsdf, sample);
        if (pathLength+1 == maxIndirectDepth) break;
        // Sample BSDF to get new path direction
        // Get random numbers for sampling new direction, \mono{bs1}, \mono{bs2}, and \mono{bcs}
        Vector wi;
        float pdf;
        BxDFType flags;
        Spectrum f = bsdf->Sample_f(wo, &wi, BSDFSample(*sample->rng),
            &pdf, BSDF_ALL, &flags);
        if (f.IsBlack() || pdf == 0.)
            break;
        specularBounce = (flags & BSDF_SPECULAR) != 0;
        pathThroughput *= f * AbsDot(wi, n) / pdf;
        ray = RayDifferential(p, wi, ray, isect.rayEpsilon);
        // Possibly terminate the path
        if (pathLength > 2) {
            float rrProb = min(1.f, pathThroughput.y());
            if (sample->rng->RandomFloat() > rrProb)
                break;
            pathThroughput /= rrProb;
        }
    }
    return L;
}

Spectrum DirectLightingIntegrator::Li(const Scene *scene,
                                      const Renderer *renderer, const RayDifferential &ray,
                                      const Intersection &isect, const Sample *sample, RNG &rng, MemoryArena &arena) const {
    Spectrum L(0.f);
    // Evaluate BSDF at hit point
    BSDF *bsdf = isect.GetBSDF(ray, arena);
    Vector wo = -ray.d;
    const Point &p = bsdf->dgShading.p;
    const Normal &n = bsdf->dgShading.nn;
    // Compute emitted light if ray hit an area light source
    L += isect.Le(wo);

    // Compute direct lighting for _DirectLightingIntegrator_ integrator
    if (scene->lights.size() > 0) {
        // Apply direct lighting strategy
        switch (strategy) {
        case SAMPLE_ALL_UNIFORM:
            L += UniformSampleAllLights(scene, renderer, arena, p, n, wo,
                                        isect.rayEpsilon, ray.time, bsdf, sample, rng,
                                        lightSampleOffsets, bsdfSampleOffsets);
            break;
        case SAMPLE_ONE_UNIFORM:
            L += UniformSampleOneLight(scene, renderer, arena, p, n, wo,
                                       isect.rayEpsilon, ray.time, bsdf, sample, rng,
                                       lightNumOffset, lightSampleOffsets, bsdfSampleOffsets);
            break;
        }
    }
    if (ray.depth + 1 < maxDepth) {
        Vector wi;
        // Trace rays for specular reflection and refraction
        Spectrum reflection(0.f);
        Spectrum refraction(0.f);
        reflection = SpecularReflect(ray, bsdf, rng, isect, renderer, scene,
                                     sample,arena);
        refraction = SpecularTransmit(ray, bsdf, rng, isect, renderer, scene,
                                      sample, arena);

        L += reflection;
        L += refraction;
    }
    //printf("Size %i \n", NaiadFoam::cur->FoamPlane().size());

    return L*bsdf->dgShading.mult;
}

int BidirIntegrator::generatePath(const Scene *scene, const Ray &r,
		const Sample *sample, const int *bsdfOffset,
		const int *bsdfCompOffset,
		BidirVertex *vertices, int maxVerts) const {
	int nVerts = 0;
	RayDifferential ray(r.o, r.d);
	while (nVerts < maxVerts) {
		// Find next vertex in path and initialize _vertices_
		Intersection isect;
		if (!scene->Intersect(ray, &isect))
			break;
		BidirVertex &v = vertices[nVerts];
		v.bsdf = isect.GetBSDF(ray); // do before Ns is set!
		v.p = isect.dg.p;
		v.ng = isect.dg.nn;
		v.ns = v.bsdf->dgShading.nn;
		v.wi = -ray.d;
		++nVerts;
		// Possibly terminate bidirectional path sampling
		if (nVerts > 2) {
			float rrProb = .2f;
			if (RandomFloat() > rrProb)
				break;
			v.rrWeight = 1.f / rrProb;
		}
		// Initialize _ray_ for next segment of path
		float u1 = sample->twoD[bsdfOffset[nVerts-1]][0];
		float u2 = sample->twoD[bsdfOffset[nVerts-1]][1];
		float u3 = sample->oneD[bsdfCompOffset[nVerts-1]][0];
		Spectrum fr = v.bsdf->Sample_f(v.wi, &v.wo, u1, u2, u3,
			 &v.bsdfWeight, BSDF_ALL, &v.flags);
		if (fr.Black() && v.bsdfWeight == 0.f)
			break;
		ray = RayDifferential(v.p, v.wo);
	}
	// Initialize additional values in _vertices_
	for (int i = 0; i < nVerts-1; ++i)
		vertices[i].dAWeight = vertices[i].bsdfWeight *
			AbsDot(-vertices[i].wo, vertices[i+1].ng) /
			DistanceSquared(vertices[i].p, vertices[i+1].p);
	return nVerts;
}

// WhittedIntegrator Method Definitions
Spectrum WhittedIntegrator::Li(const Scene *scene,
        const Renderer *renderer, const RayDifferential &ray,
        const Intersection &isect, const Sample *sample,
        MemoryArena &arena) const {
    Spectrum L(0.);
    // Compute emitted and reflected light at ray intersection point

    // Evaluate BSDF at hit point
    BSDF *bsdf = isect.GetBSDF(ray, arena);

    // Initialize common variables for Whitted integrator
    const Point &p = bsdf->dgShading.p;
    const Normal &n = bsdf->dgShading.nn;
    Vector wo = -ray.d;

    // Compute emitted light if ray hit an area light source
    L += isect.Le(wo);

    // Add contribution of each light source
    Vector wi;
    for (u_int i = 0; i < scene->lights.size(); ++i) {
        VisibilityTester visibility;
        float pdf;
        Spectrum Li = scene->lights[i]->Sample_L(p, isect.rayEpsilon,
            LightSample(*sample->rng), sample->time, &wi, &pdf, &visibility);
        if (Li.IsBlack() || pdf == 0.f) continue;
        Li /= pdf;
        Spectrum f = bsdf->f(wo, wi);
        if (!f.IsBlack() && visibility.Unoccluded(scene))
            L += f * Li * AbsDot(wi, n) *
                visibility.Transmittance(scene, renderer,
                                         sample, NULL, arena);
    }
    if (ray.depth + 1 < maxDepth) {
        // Trace rays for specular reflection and refraction
        L += SpecularReflect(ray, bsdf, *sample->rng, isect, renderer,
                             scene, sample, arena);
        L += SpecularTransmit(ray, bsdf, *sample->rng, isect, renderer,
                              scene, sample, arena);
    }
    return L;
}

// AmbientOcclusionIntegrator Method Definitions
Spectrum AmbientOcclusionIntegrator::Li(const Scene *scene, const Renderer *renderer,
        const RayDifferential &ray, const Intersection &isect,
        const Sample *sample, RNG &rng, MemoryArena &arena, int wavelength) const {

    BSDF *bsdf = isect.GetBSDF(ray, arena, wavelength);
    const Point &p = bsdf->dgShading.p;
    Normal n = Faceforward(isect.dg.nn, -ray.d);

    uint32_t scramble[2] = { rng.RandomUInt(), rng.RandomUInt() };
    float u[2];
    int nClear = 0;
    for (int i = 0; i < nSamples; ++i) {
        Sample02(i, scramble, u);
        Vector w = UniformSampleSphere(u[0], u[1]);
        if (Dot(w, n) < 0.) w = -w;
        Ray r(p, w, .01f, maxDist);
        if (!scene->IntersectP(r)) ++nClear;
    }
    return Spectrum(float(nClear) / float(nSamples));
}

Spectrum IrradianceCacheIntegrator::Li(const Scene *scene,
        const Renderer *renderer, const RayDifferential &ray, const Intersection &isect,
        const Sample *sample, RNG &rng, MemoryArena &arena) const {
    Spectrum L(0.);
    // Evaluate BSDF at hit point
    BSDF *bsdf = isect.GetBSDF(ray, arena);
    Vector wo = -ray.d;
    const Point &p = bsdf->dgShading.p;
    const Normal &n = bsdf->dgShading.nn;
    L += isect.Le(wo);
    // Compute direct lighting for irradiance cache
    L += UniformSampleAllLights(scene, renderer, arena, p, n, wo,
             isect.rayEpsilon, ray.time, bsdf, sample, rng,
             lightSampleOffsets, bsdfSampleOffsets);

    // Compute indirect lighting for irradiance cache
    if (ray.depth + 1 < maxSpecularDepth) {
        Vector wi;
        // Trace rays for specular reflection and refraction
        L += SpecularReflect(ray, bsdf, rng, isect, renderer, scene, sample,
                             arena);
        L += SpecularTransmit(ray, bsdf, rng, isect, renderer, scene, sample,
                              arena);
    }

    // Estimate indirect lighting with irradiance cache
    Normal ng = isect.dg.nn;
    ng = Faceforward(ng, wo);

    // Compute pixel spacing in world space at intersection point
    float pixelSpacing = sqrtf(Cross(isect.dg.dpdx, isect.dg.dpdy).Length());
    BxDFType flags = BxDFType(BSDF_REFLECTION | BSDF_DIFFUSE | BSDF_GLOSSY);
    L += indirectLo(p, ng, pixelSpacing, wo, isect.rayEpsilon,
                    bsdf, flags, rng, scene, renderer, arena);
    flags = BxDFType(BSDF_TRANSMISSION | BSDF_DIFFUSE | BSDF_GLOSSY);
    L += indirectLo(p, -ng, pixelSpacing, wo, isect.rayEpsilon,
                    bsdf, flags, rng, scene, renderer, arena);
    return L;
}

Spectrum IrradianceCache::Li(const Scene *scene, const RayDifferential &ray,
		const Sample *sample, float *alpha) const {
	Intersection isect;
	Spectrum L(0.);
	if (scene->Intersect(ray, &isect)) {
		if (alpha) *alpha = 1.;
		// Evaluate BSDF at hit point
		BSDF *bsdf = isect.GetBSDF(ray);
		Vector wo = -ray.d;
		const Point &p = bsdf->dgShading.p;
		const Normal &n = bsdf->dgShading.nn;
		// Compute direct lighting for irradiance cache
		L += isect.Le(wo);
		L += UniformSampleAllLights(scene, p, n, wo, bsdf, sample,
			lightSampleOffset, bsdfSampleOffset,
			bsdfComponentOffset);
		// Compute indirect lighting for irradiance cache
		if (specularDepth++ < maxSpecularDepth) {
			Vector wi;
			// Trace rays for specular reflection and refraction
			Spectrum f = bsdf->Sample_f(wo, &wi,
				BxDFType(BSDF_REFLECTION | BSDF_SPECULAR));
			if (!f.Black()) {
				// Compute ray differential _rd_ for specular reflection
				RayDifferential rd(p, wi);
				rd.hasDifferentials = true;
				rd.rx.o = p + isect.dg.dpdx;
				rd.ry.o = p + isect.dg.dpdy;
				// Compute differential reflected directions
				Normal dndx = bsdf->dgShading.dndu * bsdf->dgShading.dudx +
					bsdf->dgShading.dndv * bsdf->dgShading.dvdx;
				Normal dndy = bsdf->dgShading.dndu * bsdf->dgShading.dudy +
					bsdf->dgShading.dndv * bsdf->dgShading.dvdy;
				Vector dwodx = -ray.rx.d - wo, dwody = -ray.ry.d - wo;
				float dDNdx = Dot(dwodx, n) + Dot(wo, dndx);
				float dDNdy = Dot(dwody, n) + Dot(wo, dndy);
				rd.rx.d = wi -
				          dwodx + 2 * Vector(Dot(wo, n) * dndx +
						  dDNdx * n);
				rd.ry.d = wi -
				          dwody + 2 * Vector(Dot(wo, n) * dndy +
						  dDNdy * n);
				L += scene->Li(rd, sample) * f * AbsDot(wi, n);
			}
			f = bsdf->Sample_f(wo, &wi,
				BxDFType(BSDF_TRANSMISSION | BSDF_SPECULAR));
			if (!f.Black()) {
				// Compute ray differential _rd_ for specular transmission
				RayDifferential rd(p, wi);
				rd.hasDifferentials = true;
				rd.rx.o = p + isect.dg.dpdx;
				rd.ry.o = p + isect.dg.dpdy;
				
				float eta = bsdf->eta;
				Vector w = -wo;
				if (Dot(wo, n) < 0) eta = 1.f / eta;
				
				Normal dndx = bsdf->dgShading.dndu * bsdf->dgShading.dudx + bsdf->dgShading.dndv * bsdf->dgShading.dvdx;
				Normal dndy = bsdf->dgShading.dndu * bsdf->dgShading.dudy + bsdf->dgShading.dndv * bsdf->dgShading.dvdy;
				
				Vector dwodx = -ray.rx.d - wo, dwody = -ray.ry.d - wo;
				float dDNdx = Dot(dwodx, n) + Dot(wo, dndx);
				float dDNdy = Dot(dwody, n) + Dot(wo, dndy);
				
				float mu = eta * Dot(w, n) - Dot(wi, n);
				float dmudx = (eta - (eta*eta*Dot(w,n))/Dot(wi, n)) * dDNdx;
				float dmudy = (eta - (eta*eta*Dot(w,n))/Dot(wi, n)) * dDNdy;
				
				rd.rx.d = wi + eta * dwodx - Vector(mu * dndx + dmudx * n);
				rd.ry.d = wi + eta * dwody - Vector(mu * dndy + dmudy * n);
				L += scene->Li(rd, sample) * f * AbsDot(wi, n);
			}
		}
		--specularDepth;
		// Estimate indirect lighting with irradiance cache
		Normal ng = isect.dg.nn;
		if (Dot(wo, ng) < 0.f) ng = -ng;
		BxDFType flags = BxDFType(BSDF_REFLECTION |
		                          BSDF_DIFFUSE |
								  BSDF_GLOSSY);
		L += IndirectLo(p, ng, wo, bsdf, flags, sample, scene);
		flags = BxDFType(BSDF_TRANSMISSION |
		                 BSDF_DIFFUSE |
						 BSDF_GLOSSY);
		L += IndirectLo(p, -ng, wo, bsdf, flags, sample, scene);
	}
	else {
		// Handle ray with no intersection
		if (alpha) *alpha = 0.;
		for (u_int i = 0; i < scene->lights.size(); ++i)
			L += scene->lights[i]->Le(ray);
		if (alpha && !L.Black()) *alpha = 1.;
		return L;
	}
	return L;
}

Spectrum PhotonIntegrator::Li(const Scene *scene, const Renderer *renderer,
        const RayDifferential &ray, const Intersection &isect,
        const Sample *sample, RNG &rng, MemoryArena &arena) const {
    Spectrum L(0.);
    Vector wo = -ray.d;
    // Compute emitted light if ray hit an area light source
    L += isect.Le(wo);

    // Evaluate BSDF at hit pbrt::Point
    BSDF *bsdf = isect.GetBSDF(ray, arena);
    const pbrt::Point &p = bsdf->dgShading.p;
    const Normal &n = bsdf->dgShading.nn;
    L += UniformSampleAllLights(scene, renderer, arena, p, n,
        wo, isect.rayEpsilon, ray.time, bsdf, sample, rng,
        lightSampleOffsets, bsdfSampleOffsets);
    // Compute caustic lighting for photon map integrator
    ClosePhoton *lookupBuf = arena.Alloc<ClosePhoton>(nLookup);
    L += LPhoton(causticMap, nCausticPaths, nLookup, lookupBuf, bsdf,
                 rng, isect, wo, maxDistSquared);

    // Compute indirect lighting for photon map integrator
    if (finalGather && indirectMap != NULL) {
    #if 1
        // Do one-bounce final gather for photon map
        BxDFType nonSpecular = BxDFType(BSDF_REFLECTION |
            BSDF_TRANSMISSION | BSDF_DIFFUSE | BSDF_GLOSSY);
        if (bsdf->NumComponents(nonSpecular) > 0) {
            // Find indirect photons around point for importance sampling
            const uint32_t nIndirSamplePhotons = 50;
            PhotonProcess proc(nIndirSamplePhotons,
                               arena.Alloc<ClosePhoton>(nIndirSamplePhotons));
            float searchDist2 = maxDistSquared;
            while (proc.nFound < nIndirSamplePhotons) {
                float md2 = searchDist2;
                proc.nFound = 0;
                indirectMap->Lookup(p, proc, md2);
                searchDist2 *= 2.f;
            }

            // Copy photon directions to local array
            Vector *photonDirs = arena.Alloc<Vector>(nIndirSamplePhotons);
            for (uint32_t i = 0; i < nIndirSamplePhotons; ++i)
                photonDirs[i] = proc.photons[i].photon->wi;

            // Use BSDF to do final gathering
            Spectrum Li = 0.;
            for (int i = 0; i < gatherSamples; ++i) {
                // Sample random direction from BSDF for final gather ray
                Vector wi;
                float pdf;
                BSDFSample bsdfSample(sample, bsdfGatherSampleOffsets, i);
                Spectrum fr = bsdf->Sample_f(wo, &wi, bsdfSample,
                                             &pdf, BxDFType(BSDF_ALL & ~BSDF_SPECULAR));
                if (fr.IsBlack() || pdf == 0.f) continue;
                Assert(pdf >= 0.f);

                // Trace BSDF final gather ray and accumulate radiance
                RayDifferential bounceRay(p, wi, ray, isect.rayEpsilon);
                Intersection gatherIsect;
                if (scene->Intersect(bounceRay, &gatherIsect)) {
                    // Compute exitant radiance _Lindir_ using radiance photons
                    Spectrum Lindir = 0.f;
                    Normal nGather = gatherIsect.dg.nn;
                    nGather = Faceforward(nGather, -bounceRay.d);
                    RadiancePhotonProcess proc(nGather);
                    float md2 = INFINITY;
                    radianceMap->Lookup(gatherIsect.dg.p, proc, md2);
                    if (proc.photon != NULL)
                        Lindir = proc.photon->Lo;
                    Lindir *= renderer->Transmittance(scene, bounceRay, NULL, rng, arena);

                    // Compute MIS weight for BSDF-sampled gather ray

                    // Compute PDF for photon-sampling of direction _wi_
                    float photonPdf = 0.f;
                    float conePdf = UniformConePdf(cosGatherAngle);
                    for (uint32_t j = 0; j < nIndirSamplePhotons; ++j)
                        if (Dot(photonDirs[j], wi) > .999f * cosGatherAngle)
                            photonPdf += conePdf;
                    photonPdf /= nIndirSamplePhotons;
                    float wt = PowerHeuristic(gatherSamples, pdf, gatherSamples, photonPdf);
                    Li += fr * Lindir * (AbsDot(wi, n) * wt / pdf);
                }
            }
            L += Li / gatherSamples;

            // Use nearby photons to do final gathering
            Li = 0.;
            for (int i = 0; i < gatherSamples; ++i) {
                // Sample random direction using photons for final gather ray
                BSDFSample gatherSample(sample, indirGatherSampleOffsets, i);
                int photonNum = min((int)nIndirSamplePhotons - 1,
                    Floor2Int(gatherSample.uComponent * nIndirSamplePhotons));

                // Sample gather ray direction from _photonNum_
                Vector vx, vy;
                CoordinateSystem(photonDirs[photonNum], &vx, &vy);
                Vector wi = UniformSampleCone(gatherSample.uDir[0], gatherSample.uDir[1],
                                              cosGatherAngle, vx, vy, photonDirs[photonNum]);

//.........这里部分代码省略.........

void PhotonShootingTask::Run() {
    // Declare local variables for _PhotonShootingTask_
    MemoryArena arena;
    RNG rng(31 * taskNum);
    vector<Photon> localDirectPhotons, localIndirectPhotons, localCausticPhotons;
    vector<RadiancePhoton> localRadiancePhotons;
    uint32_t totalPaths = 0;
    bool causticDone = (integrator->nCausticPhotonsWanted == 0);
    bool indirectDone = (integrator->nIndirectPhotonsWanted == 0);
    PermutedHalton halton(6, rng);
    vector<Spectrum> localRpReflectances, localRpTransmittances;
    while (true) {
        // Follow photon paths for a block of samples
        const uint32_t blockSize = 4096;
        for (uint32_t i = 0; i < blockSize; ++i) {
            float u[6];
            halton.Sample(++totalPaths, u);
            // Choose light to shoot photon from
            float lightPdf;
            int lightNum = lightDistribution->SampleDiscrete(u[0], &lightPdf);
            const Light *light = scene->lights[lightNum];

            // Generate _photonRay_ from light source and initialize _alpha_
            RayDifferential photonRay;
            float pdf;
            LightSample ls(u[1], u[2], u[3]);
            Normal Nl;
            Spectrum Le = light->Sample_L(scene, ls, u[4], u[5],
                                          time, &photonRay, &Nl, &pdf);
            if (pdf == 0.f || Le.IsBlack()) continue;
            Spectrum alpha = (AbsDot(Nl, photonRay.d) * Le) / (pdf * lightPdf);
            if (!alpha.IsBlack()) {
                // Follow photon path through scene and record intersections
                PBRT_PHOTON_MAP_STARTED_RAY_PATH(&photonRay, &alpha);
                bool specularPath = true;
                Intersection photonIsect;
                int nIntersections = 0;
                while (scene->Intersect(photonRay, &photonIsect)) {
                    ++nIntersections;
                    // Handle photon/surface intersection
                    alpha *= renderer->Transmittance(scene, photonRay, NULL, rng, arena);
                    BSDF *photonBSDF = photonIsect.GetBSDF(photonRay, arena);
                    BxDFType specularType = BxDFType(BSDF_REFLECTION |
                                            BSDF_TRANSMISSION | BSDF_SPECULAR);
                    bool hasNonSpecular = (photonBSDF->NumComponents() >
                                           photonBSDF->NumComponents(specularType));
                    Vector wo = -photonRay.d;
                    if (hasNonSpecular) {
                        // Deposit photon at surface
                        Photon photon(photonIsect.dg.p, alpha, wo);
                        bool depositedPhoton = false;
                        if (specularPath && nIntersections > 1) {
                            if (!causticDone) {
                                PBRT_PHOTON_MAP_DEPOSITED_CAUSTIC_PHOTON(&photonIsect.dg, &alpha, &wo);
                                depositedPhoton = true;
                                localCausticPhotons.push_back(photon);
                            }
                        }
                        else {
                            // Deposit either direct or indirect photon
                            // stop depositing direct photons once indirectDone is true; don't
                            // want to waste memory storing too many if we're going a long time
                            // trying to get enough caustic photons desposited.
                            if (nIntersections == 1 && !indirectDone && integrator->finalGather) {
                                PBRT_PHOTON_MAP_DEPOSITED_DIRECT_PHOTON(&photonIsect.dg, &alpha, &wo);
                                depositedPhoton = true;
                                localDirectPhotons.push_back(photon);
                            }
                            else if (nIntersections > 1 && !indirectDone) {
                                PBRT_PHOTON_MAP_DEPOSITED_INDIRECT_PHOTON(&photonIsect.dg, &alpha, &wo);
                                depositedPhoton = true;
                                localIndirectPhotons.push_back(photon);
                            }
                        }

                        // Possibly create radiance photon at photon intersection point
                        if (depositedPhoton && integrator->finalGather &&
                                rng.RandomFloat() < .125f) {
                            Normal n = photonIsect.dg.nn;
                            n = Faceforward(n, -photonRay.d);
                            localRadiancePhotons.push_back(RadiancePhoton(photonIsect.dg.p, n));
                            Spectrum rho_r = photonBSDF->rho(rng, BSDF_ALL_REFLECTION);
                            localRpReflectances.push_back(rho_r);
                            Spectrum rho_t = photonBSDF->rho(rng, BSDF_ALL_TRANSMISSION);
                            localRpTransmittances.push_back(rho_t);
                        }
                    }
                    if (nIntersections >= integrator->maxPhotonDepth) break;

                    // Sample new photon ray direction
                    Vector wi;
                    float pdf;
                    BxDFType flags;
                    Spectrum fr = photonBSDF->Sample_f(wo, &wi, BSDFSample(rng),
                                                       &pdf, BSDF_ALL, &flags);
                    if (fr.IsBlack() || pdf == 0.f) break;
                    Spectrum anew = alpha * fr *
                        AbsDot(wi, photonBSDF->dgShading.nn) / pdf;

                    // Possibly terminate photon path with Russian roulette
//.........这里部分代码省略.........

void PhotonIntegrator::Preprocess(const Scene *scene) {
	if (scene->lights.size() == 0) return;
	ProgressReporter progress(nCausticPhotons+nDirectPhotons+ // NOBOOK
		nIndirectPhotons, "Shooting photons"); // NOBOOK
	vector<Photon> causticPhotons;
	vector<Photon> directPhotons;
	vector<Photon> indirectPhotons;
	causticPhotons.reserve(nCausticPhotons); // NOBOOK
	directPhotons.reserve(nDirectPhotons); // NOBOOK
	indirectPhotons.reserve(nIndirectPhotons); // NOBOOK
	// Initialize photon shooting statistics
	static StatsCounter nshot("Photon Map",
		"Number of photons shot from lights");
	bool causticDone = (nCausticPhotons == 0);
	bool directDone = (nDirectPhotons == 0);
	bool indirectDone = (nIndirectPhotons == 0);
	while (!causticDone || !directDone || !indirectDone) {
		++nshot;
		// Give up if we're not storing enough photons
		if (nshot > 500000 &&
			(unsuccessful(nCausticPhotons,
			              causticPhotons.size(),
						  nshot) ||
			 unsuccessful(nDirectPhotons,
			              directPhotons.size(),
						  nshot) ||
			 unsuccessful(nIndirectPhotons,
			              indirectPhotons.size(),
						  nshot))) {
			Error("Unable to store enough photons.  Giving up.\n");
			return;
		}
		// Trace a photon path and store contribution
		// Choose 4D sample values for photon
		float u[4];
		u[0] = (float)RadicalInverse((int)nshot+1, 2);
		u[1] = (float)RadicalInverse((int)nshot+1, 3);
		u[2] = (float)RadicalInverse((int)nshot+1, 5);
		u[3] = (float)RadicalInverse((int)nshot+1, 7);
		// Choose light to shoot photon from
		int nLights = int(scene->lights.size());
		int lightNum =
			min(Floor2Int(nLights * (float)RadicalInverse((int)nshot+1, 11)),
			nLights-1);
		Light *light = scene->lights[lightNum];
		float lightPdf = 1.f / nLights;
		// Generate _photonRay_ from light source and initialize _alpha_
		RayDifferential photonRay;
		float pdf;
		Spectrum alpha =
			light->Sample_L(scene, u[0], u[1], u[2], u[3],
				&photonRay, &pdf);
		if (pdf == 0.f || alpha.Black()) continue;
		alpha /= pdf * lightPdf;
		if (!alpha.Black()) {
			// Follow photon path through scene and record intersections
			bool specularPath = false;
			Intersection photonIsect;
			int nIntersections = 0;
			while (scene->Intersect(photonRay, &photonIsect)) {
				++nIntersections;
				// Handle photon/surface intersection
				alpha *= scene->Transmittance(photonRay);
				Vector wo = -photonRay.d;
				BSDF *photonBSDF = photonIsect.GetBSDF(photonRay);
				BxDFType specularType = BxDFType(BSDF_REFLECTION |
					BSDF_TRANSMISSION | BSDF_SPECULAR);
				bool hasNonSpecular = (photonBSDF->NumComponents() >
					photonBSDF->NumComponents(specularType));
				if (hasNonSpecular) {
					// Deposit photon at surface
					Photon photon(photonIsect.dg.p, alpha, wo);
					if (nIntersections == 1) {
						// Process direct lighting photon intersection
						if (!directDone) {
							directPhotons.push_back(photon);
							if (directPhotons.size() == nDirectPhotons) {
								directDone = true;
								nDirectPaths = (int)nshot;
								directMap =
									new KdTree<Photon,
											   PhotonProcess>(directPhotons);
							}
							progress.Update(); // NOBOOK
						}
					}
					else if (specularPath) {
						// Process caustic photon intersection
						if (!causticDone) {
							causticPhotons.push_back(photon);
							if (causticPhotons.size() == nCausticPhotons) {
								causticDone = true;
								nCausticPaths = (int)nshot;
								causticMap =
									new KdTree<Photon,
										       PhotonProcess>(causticPhotons);
							}
							progress.Update();
						}
					}
//.........这里部分代码省略.........

// WhittedIntegrator Method Definitions
Spectrum WhittedIntegrator::Li(const Scene *scene,
		const RayDifferential &ray, const Sample *sample,
		float *alpha) const {
	Intersection isect;
	Spectrum L(0.);
	bool hitSomething;
	// Search for ray-primitive intersection
	hitSomething = scene->Intersect(ray, &isect);
	if (!hitSomething) {
		// Handle ray with no intersection
		if (alpha) *alpha = 0.;
		for (u_int i = 0; i < scene->lights.size(); ++i)
			L += scene->lights[i]->Le(ray);
		if (alpha && !L.Black()) *alpha = 1.;
		return L;
	}
	else {
		// Initialize _alpha_ for ray hit
		if (alpha) *alpha = 1.;
		// Compute emitted and reflected light at ray intersection point
		// Evaluate BSDF at hit point
		BSDF *bsdf = isect.GetBSDF(ray);
		// Initialize common variables for Whitted integrator
		const Point &p = bsdf->dgShading.p;
		const Normal &n = bsdf->dgShading.nn;
		Vector wo = -ray.d;
		// Compute emitted light if ray hit an area light source
		L += isect.Le(wo);
		// Add contribution of each light source
		Vector wi;
		for (u_int i = 0; i < scene->lights.size(); ++i) {
			VisibilityTester visibility;
			Spectrum Li = scene->lights[i]->Sample_L(p, &wi, &visibility);
			if (Li.Black()) continue;
			Spectrum f = bsdf->f(wo, wi);
			if (!f.Black() && visibility.Unoccluded(scene))
				L += f * Li * AbsDot(wi, n) * visibility.Transmittance(scene);
		}
		if (rayDepth++ < maxDepth) {
			// Trace rays for specular reflection and refraction
			Spectrum f = bsdf->Sample_f(wo, &wi,
				BxDFType(BSDF_REFLECTION | BSDF_SPECULAR));
			if (!f.Black() && AbsDot(wi, n) > 0.f) {
				// Compute ray differential _rd_ for specular reflection
				RayDifferential rd(p, wi);
				rd.hasDifferentials = true;
				rd.rx.o = p + isect.dg.dpdx;
				rd.ry.o = p + isect.dg.dpdy;
				// Compute differential reflected directions
				Normal dndx = bsdf->dgShading.dndu * bsdf->dgShading.dudx +
					bsdf->dgShading.dndv * bsdf->dgShading.dvdx;
				Normal dndy = bsdf->dgShading.dndu * bsdf->dgShading.dudy +
					bsdf->dgShading.dndv * bsdf->dgShading.dvdy;
				Vector dwodx = -ray.rx.d - wo, dwody = -ray.ry.d - wo;
				float dDNdx = Dot(dwodx, n) + Dot(wo, dndx);
				float dDNdy = Dot(dwody, n) + Dot(wo, dndy);
				rd.rx.d = wi -
				          dwodx + 2 * Vector(Dot(wo, n) * dndx +
						  dDNdx * n);
				rd.ry.d = wi -
				          dwody + 2 * Vector(Dot(wo, n) * dndy +
						  dDNdy * n);
				L += scene->Li(rd, sample) * f * AbsDot(wi, n);
			}
			f = bsdf->Sample_f(wo, &wi,
				BxDFType(BSDF_TRANSMISSION | BSDF_SPECULAR));
			if (!f.Black() && AbsDot(wi, n) > 0.f) {
				// Compute ray differential _rd_ for specular transmission
				RayDifferential rd(p, wi);
				rd.hasDifferentials = true;
				rd.rx.o = p + isect.dg.dpdx;
				rd.ry.o = p + isect.dg.dpdy;
				
				float eta = bsdf->eta;
				Vector w = -wo;
				if (Dot(wo, n) < 0) eta = 1.f / eta;
				
				Normal dndx = bsdf->dgShading.dndu * bsdf->dgShading.dudx + bsdf->dgShading.dndv * bsdf->dgShading.dvdx;
				Normal dndy = bsdf->dgShading.dndu * bsdf->dgShading.dudy + bsdf->dgShading.dndv * bsdf->dgShading.dvdy;
				
				Vector dwodx = -ray.rx.d - wo, dwody = -ray.ry.d - wo;
				float dDNdx = Dot(dwodx, n) + Dot(wo, dndx);
				float dDNdy = Dot(dwody, n) + Dot(wo, dndy);
				
				float mu = eta * Dot(w, n) - Dot(wi, n);
				float dmudx = (eta - (eta*eta*Dot(w,n))/Dot(wi, n)) * dDNdx;
				float dmudy = (eta - (eta*eta*Dot(w,n))/Dot(wi, n)) * dDNdy;
				
				rd.rx.d = wi + eta * dwodx - Vector(mu * dndx + dmudx * n);
				rd.ry.d = wi + eta * dwody - Vector(mu * dndy + dmudy * n);
				L += scene->Li(rd, sample) * f * AbsDot(wi, n);
			}
		}
		--rayDepth;
	}
	return L;
}

Spectrum PhotonIntegrator::Li(const Scene *scene,
		const RayDifferential &ray, const Sample *sample,
		float *alpha) const {
	// Compute reflected radiance with photon map
	Spectrum L(0.);
	Intersection isect;
	if (scene->Intersect(ray, &isect)) {
		if (alpha) *alpha = 1.;
		Vector wo = -ray.d;
		// Compute emitted light if ray hit an area light source
		L += isect.Le(wo);
		// Evaluate BSDF at hit point
		BSDF *bsdf = isect.GetBSDF(ray);
		const Point &p = bsdf->dgShading.p;
		const Normal &n = bsdf->dgShading.nn;
		// Compute direct lighting for photon map integrator
		if (directWithPhotons)
			L += LPhoton(directMap, nDirectPaths, nLookup,
				bsdf, isect, wo, maxDistSquared);
		else
			L += UniformSampleAllLights(scene, p, n,
				wo, bsdf, sample,
				lightSampleOffset, bsdfSampleOffset,
				bsdfComponentOffset);
		
		// Compute indirect lighting for photon map integrator
		L += LPhoton(causticMap, nCausticPaths, nLookup, bsdf,
			isect, wo, maxDistSquared);
		if (finalGather) {
			// Do one-bounce final gather for photon map
			Spectrum Li(0.);
			for (int i = 0; i < gatherSamples; ++i) {
				// Sample random direction for final gather ray
				Vector wi;
				float u1 = sample->twoD[gatherSampleOffset][2*i];
				float u2 = sample->twoD[gatherSampleOffset][2*i+1];
				float u3 = sample->oneD[gatherComponentOffset][i];
				float pdf;
				Spectrum fr = bsdf->Sample_f(wo, &wi, u1, u2, u3,
					&pdf, BxDFType(BSDF_ALL & (~BSDF_SPECULAR)));
				if (fr.Black() || pdf == 0.f) continue;
				RayDifferential bounceRay(p, wi);
				static StatsCounter gatherRays("Photon Map", // NOBOOK
					"Final gather rays traced"); // NOBOOK
				++gatherRays; // NOBOOK
				Intersection gatherIsect;
				if (scene->Intersect(bounceRay, &gatherIsect)) {
					// Compute exitant radiance at final gather intersection
					BSDF *gatherBSDF = gatherIsect.GetBSDF(bounceRay);
					Vector bounceWo = -bounceRay.d;
					Spectrum Lindir =
						LPhoton(directMap, nDirectPaths, nLookup,
							gatherBSDF, gatherIsect, bounceWo, maxDistSquared) +
						LPhoton(indirectMap, nIndirectPaths, nLookup,
							gatherBSDF, gatherIsect, bounceWo, maxDistSquared) +
						LPhoton(causticMap, nCausticPaths, nLookup,
							gatherBSDF, gatherIsect, bounceWo, maxDistSquared);
					Lindir *= scene->Transmittance(bounceRay);
					Li += fr * Lindir * AbsDot(wi, n) / pdf;
				}
			}
			L += Li / float(gatherSamples);
		}
		else
			L += LPhoton(indirectMap, nIndirectPaths, nLookup,
				bsdf, isect, wo, maxDistSquared);
		if (specularDepth++ < maxSpecularDepth) {
			Vector wi;
			// Trace rays for specular reflection and refraction
			Spectrum f = bsdf->Sample_f(wo, &wi,
				BxDFType(BSDF_REFLECTION | BSDF_SPECULAR));
			if (!f.Black()) {
				// Compute ray differential _rd_ for specular reflection
				RayDifferential rd(p, wi);
				rd.hasDifferentials = true;
				rd.rx.o = p + isect.dg.dpdx;
				rd.ry.o = p + isect.dg.dpdy;
				// Compute differential reflected directions
				Normal dndx = bsdf->dgShading.dndu * bsdf->dgShading.dudx +
					bsdf->dgShading.dndv * bsdf->dgShading.dvdx;
				Normal dndy = bsdf->dgShading.dndu * bsdf->dgShading.dudy +
					bsdf->dgShading.dndv * bsdf->dgShading.dvdy;
				Vector dwodx = -ray.rx.d - wo, dwody = -ray.ry.d - wo;
				float dDNdx = Dot(dwodx, n) + Dot(wo, dndx);
				float dDNdy = Dot(dwody, n) + Dot(wo, dndy);
				rd.rx.d = wi -
				          dwodx + 2 * Vector(Dot(wo, n) * dndx +
						  dDNdx * n);
				rd.ry.d = wi -
				          dwody + 2 * Vector(Dot(wo, n) * dndy +
						  dDNdy * n);
				L += scene->Li(rd, sample) * f * AbsDot(wi, n);
			}
			f = bsdf->Sample_f(wo, &wi,
				BxDFType(BSDF_TRANSMISSION | BSDF_SPECULAR));
			if (!f.Black()) {
				// Compute ray differential _rd_ for specular transmission
				RayDifferential rd(p, wi);
				rd.hasDifferentials = true;
				rd.rx.o = p + isect.dg.dpdx;
//.........这里部分代码省略.........