00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025 #include "pbrt.h"
00026 #include "transport.h"
00027 #include "scene.h"
00028 #include "mc.h"
00029 #include "sampling.h"
00030
00031
00032 struct VirtualLight {
00033 VirtualLight() { }
00034 VirtualLight(const Point &pp, const Normal &nn, const Spectrum &le)
00035 : p(pp), n(nn), Le(le) { }
00036 Point p;
00037 Normal n;
00038 Spectrum Le;
00039 };
00040
00041 class IGIIntegrator : public SurfaceIntegrator {
00042 public:
00043
00044 IGIIntegrator(int nl, int ns, float md, float rrt, float is);
00045 Spectrum Li(const Scene *scene, const RayDifferential &ray,
00046 const Sample *sample, float *alpha) const;
00047 void RequestSamples(Sample *sample, const Scene *scene);
00048 void Preprocess(const Scene *);
00049
00050 private:
00051
00052 u_int nLightPaths, nLightSets;
00053 vector<VirtualLight> *virtualLights;
00054 mutable int specularDepth;
00055 int maxSpecularDepth;
00056 float minDist2, rrThreshold, indirectScale;
00057 int vlSetOffset;
00058
00059 int *lightSampleOffset, lightNumOffset;
00060 int *bsdfSampleOffset, *bsdfComponentOffset;
00061 };
00062
00063
00064 IGIIntegrator::IGIIntegrator(int nl, int ns, float md,
00065 float rrt, float is) {
00066 nLightPaths = RoundUpPow2((u_int)nl);
00067 nLightSets = RoundUpPow2((u_int)ns);
00068 minDist2 = md * md;
00069 rrThreshold = rrt;
00070 indirectScale = is;
00071 maxSpecularDepth = 5;
00072 specularDepth = 0;
00073 virtualLights = new vector<VirtualLight>[nLightSets];
00074 }
00075 void IGIIntegrator::RequestSamples(Sample *sample,
00076 const Scene *scene) {
00077
00078 u_int nLights = scene->lights.size();
00079 lightSampleOffset = new int[nLights];
00080 bsdfSampleOffset = new int[nLights];
00081 bsdfComponentOffset = new int[nLights];
00082 for (u_int i = 0; i < nLights; ++i) {
00083 const Light *light = scene->lights[i];
00084 int lightSamples =
00085 scene->sampler->RoundSize(light->nSamples);
00086 lightSampleOffset[i] = sample->Add2D(lightSamples);
00087 bsdfSampleOffset[i] = sample->Add2D(lightSamples);
00088 bsdfComponentOffset[i] = sample->Add1D(lightSamples);
00089 }
00090 lightNumOffset = -1;
00091 vlSetOffset = sample->Add1D(1);
00092 }
00093 void IGIIntegrator::Preprocess(const Scene *scene) {
00094 if (scene->lights.size() == 0) return;
00095
00096 float *lightNum = new float[nLightPaths * nLightSets];
00097 float *lightSamp0 = new float[2 * nLightPaths * nLightSets];
00098 float *lightSamp1 = new float[2 * nLightPaths * nLightSets];
00099 LDShuffleScrambled1D(nLightPaths, nLightSets, lightNum);
00100 LDShuffleScrambled2D(nLightPaths, nLightSets, lightSamp0);
00101 LDShuffleScrambled2D(nLightPaths, nLightSets, lightSamp1);
00102
00103 int nLights = int(scene->lights.size());
00104 float *lightPower = (float *)alloca(nLights * sizeof(float));
00105 float *lightCDF = (float *)alloca((nLights+1) * sizeof(float));
00106 for (int i = 0; i < nLights; ++i)
00107 lightPower[i] = scene->lights[i]->Power(scene).y();
00108 float totalPower;
00109 ComputeStep1dCDF(lightPower, nLights, &totalPower, lightCDF);
00110 for (u_int s = 0; s < nLightSets; ++s) {
00111 for (u_int i = 0; i < nLightPaths; ++i) {
00112
00113 int sampOffset = s*nLightPaths + i;
00114
00115 float lightPdf;
00116 int lNum = Floor2Int(SampleStep1d(lightPower, lightCDF,
00117 totalPower, nLights, lightNum[sampOffset], &lightPdf) * nLights);
00118 fprintf(stderr, "samp %f -> num %d\n", lightNum[sampOffset], lNum);
00119 Light *light = scene->lights[lNum];
00120
00121 RayDifferential ray;
00122 float pdf;
00123 Spectrum alpha =
00124 light->Sample_L(scene, lightSamp0[2*sampOffset],
00125 lightSamp0[2*sampOffset+1],
00126 lightSamp1[2*sampOffset],
00127 lightSamp1[2*sampOffset+1],
00128 &ray, &pdf);
00129 if (pdf == 0.f || alpha.Black()) continue;
00130 alpha /= pdf * lightPdf;
00131 fprintf(stderr, "initial alpha %f, light # %d\n", alpha.y(), lNum);
00132 Intersection isect;
00133 int nIntersections = 0;
00134 while (scene->Intersect(ray, &isect) && !alpha.Black()) {
00135 ++nIntersections;
00136 alpha *= scene->Transmittance(ray);
00137 Vector wo = -ray.d;
00138 BSDF *bsdf = isect.GetBSDF(ray);
00139
00140 static StatsCounter vls("IGI Integrator", "Virtual Lights Created");
00141 ++vls;
00142 Spectrum Le = alpha * bsdf->rho(wo) / M_PI;
00143 fprintf(stderr, "\tmade light with le y %f\n", Le.y());
00144 virtualLights[s].push_back(VirtualLight(isect.dg.p, isect.dg.nn, Le));
00145
00146 Vector wi;
00147 float pdf;
00148 BxDFType flags;
00149 Spectrum fr = bsdf->Sample_f(wo, &wi, RandomFloat(),
00150 RandomFloat(), RandomFloat(),
00151 &pdf, BSDF_ALL, &flags);
00152 if (fr.Black() || pdf == 0.f)
00153 break;
00154 Spectrum anew = alpha * fr * AbsDot(wi, bsdf->dgShading.nn) / pdf;
00155 float r = anew.y() / alpha.y();
00156 fprintf(stderr, "\tr = %f\n", r);
00157 if (RandomFloat() > r)
00158 break;
00159 alpha = anew / r;
00160 fprintf(stderr, "\tnew alpha %f\n", alpha.y());
00161 ray = RayDifferential(isect.dg.p, wi);
00162 }
00163 BSDF::FreeAll();
00164 }
00165 }
00166 delete[] lightNum;
00167 delete[] lightSamp0;
00168 delete[] lightSamp1;
00169 }
00170 Spectrum IGIIntegrator::Li(const Scene *scene,
00171 const RayDifferential &ray, const Sample *sample,
00172 float *alpha) const {
00173 Spectrum L(0.);
00174 Intersection isect;
00175 if (scene->Intersect(ray, &isect)) {
00176 if (alpha) *alpha = 1.;
00177 Vector wo = -ray.d;
00178
00179 L += isect.Le(wo);
00180
00181 BSDF *bsdf = isect.GetBSDF(ray);
00182 const Point &p = bsdf->dgShading.p;
00183 const Normal &n = bsdf->dgShading.nn;
00184 L += UniformSampleAllLights(scene, p, n,
00185 wo, bsdf, sample,
00186 lightSampleOffset, bsdfSampleOffset,
00187 bsdfComponentOffset);
00188
00189 u_int lSet = min(u_int(sample->oneD[vlSetOffset][0] * nLightSets),
00190 nLightSets-1);
00191 for (u_int i = 0; i < virtualLights[lSet].size(); ++i) {
00192 const VirtualLight &vl = virtualLights[lSet][i];
00193
00194
00195 float d2 = DistanceSquared(p, vl.p);
00196
00197 float distScale = SmoothStep(.8 * minDist2, 1.2 * minDist2, d2);
00198
00199 Vector wi = Normalize(vl.p - p);
00200 Spectrum f = distScale * bsdf->f(wo, wi);
00201 if (f.Black()) continue;
00202 float G = AbsDot(wi, n) * AbsDot(wi, vl.n) / d2;
00203 Spectrum Llight = indirectScale * f * G * vl.Le /
00204 virtualLights[lSet].size();
00205 Llight *= scene->Transmittance(Ray(p, vl.p - p));
00206
00207 if (Llight.y() < rrThreshold) {
00208 float continueProbability = .1f;
00209 if (RandomFloat() > continueProbability)
00210 continue;
00211 Llight /= continueProbability;
00212 }
00213 static StatsCounter vlsr("IGI Integrator", "Shadow Rays to Virtual Lights");
00214 ++vlsr;
00215 if (!scene->IntersectP(Ray(p, vl.p - p, RAY_EPSILON,
00216 1.f - RAY_EPSILON)))
00217 L += Llight;
00218 }
00219
00220 if (specularDepth++ < maxSpecularDepth) {
00221 Vector wi;
00222
00223 Spectrum f = bsdf->Sample_f(wo, &wi,
00224 BxDFType(BSDF_REFLECTION | BSDF_SPECULAR));
00225 if (!f.Black()) {
00226
00227 RayDifferential rd(p, wi);
00228 rd.hasDifferentials = true;
00229 rd.rx.o = p + isect.dg.dpdx;
00230 rd.ry.o = p + isect.dg.dpdy;
00231
00232 Normal dndx = bsdf->dgShading.dndu * bsdf->dgShading.dudx +
00233 bsdf->dgShading.dndv * bsdf->dgShading.dvdx;
00234 Normal dndy = bsdf->dgShading.dndu * bsdf->dgShading.dudy +
00235 bsdf->dgShading.dndv * bsdf->dgShading.dvdy;
00236 Vector dwodx = -ray.rx.d - wo, dwody = -ray.ry.d - wo;
00237 float dDNdx = Dot(dwodx, n) + Dot(wo, dndx);
00238 float dDNdy = Dot(dwody, n) + Dot(wo, dndy);
00239 rd.rx.d = wi -
00240 dwodx + 2 * Vector(Dot(wo, n) * dndx +
00241 dDNdx * n);
00242 rd.ry.d = wi -
00243 dwody + 2 * Vector(Dot(wo, n) * dndy +
00244 dDNdy * n);
00245 L += scene->Li(rd, sample) * f * AbsDot(wi, n);
00246 }
00247 f = bsdf->Sample_f(wo, &wi,
00248 BxDFType(BSDF_TRANSMISSION | BSDF_SPECULAR));
00249 if (!f.Black()) {
00250
00251 RayDifferential rd(p, wi);
00252 rd.hasDifferentials = true;
00253 rd.rx.o = p + isect.dg.dpdx;
00254 rd.ry.o = p + isect.dg.dpdy;
00255
00256 float eta = bsdf->eta;
00257 Vector w = -wo;
00258 if (Dot(wo, n) < 0) eta = 1.f / eta;
00259
00260 Normal dndx = bsdf->dgShading.dndu * bsdf->dgShading.dudx + bsdf->dgShading.dndv * bsdf->dgShading.dvdx;
00261 Normal dndy = bsdf->dgShading.dndu * bsdf->dgShading.dudy + bsdf->dgShading.dndv * bsdf->dgShading.dvdy;
00262
00263 Vector dwodx = -ray.rx.d - wo, dwody = -ray.ry.d - wo;
00264 float dDNdx = Dot(dwodx, n) + Dot(wo, dndx);
00265 float dDNdy = Dot(dwody, n) + Dot(wo, dndy);
00266
00267 float mu = eta * Dot(w, n) - Dot(wi, n);
00268 float dmudx = (eta - (eta*eta*Dot(w,n))/Dot(wi, n)) * dDNdx;
00269 float dmudy = (eta - (eta*eta*Dot(w,n))/Dot(wi, n)) * dDNdy;
00270
00271 rd.rx.d = wi + eta * dwodx - Vector(mu * dndx + dmudx * n);
00272 rd.ry.d = wi + eta * dwody - Vector(mu * dndy + dmudy * n);
00273 L += scene->Li(rd, sample) * f * AbsDot(wi, n);
00274 }
00275 }
00276 --specularDepth;
00277 }
00278 else {
00279
00280 if (alpha) *alpha = 0.;
00281 for (u_int i = 0; i < scene->lights.size(); ++i)
00282 L += scene->lights[i]->Le(ray);
00283 if (alpha && !L.Black()) *alpha = 1.;
00284 return L;
00285 }
00286 return L;
00287 }
00288 extern "C" DLLEXPORT SurfaceIntegrator *CreateSurfaceIntegrator(const ParamSet ¶ms)
00289 {
00290 int nLightPaths = params.FindOneInt("nlights", 64);
00291 int nLightSets = params.FindOneInt("nsets", 4);
00292 float minDist = params.FindOneFloat("mindist", .1f);
00293 float rrThresh = params.FindOneFloat("rrthreshold", .05f);
00294 float is = params.FindOneFloat("indirectscale", 1.f);
00295 return new IGIIntegrator(nLightPaths, nLightSets, minDist, rrThresh, is);
00296 }