00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024 #if defined(_MSC_VER)
00025 #pragma once
00026 #endif
00027
00028 #ifndef PBRT_CORE_MIPMAP_H
00029 #define PBRT_CORE_MIPMAP_H
00030
00031
00032 #include "pbrt.h"
00033 #include "spectrum.h"
00034 #include "texture.h"
00035
00036
00037 typedef enum {
00038 TEXTURE_REPEAT,
00039 TEXTURE_BLACK,
00040 TEXTURE_CLAMP
00041 } ImageWrap;
00042 template <typename T> class MIPMap {
00043 public:
00044
00045 MIPMap() { pyramid = NULL; width = height = nLevels = 0; }
00046 MIPMap(uint32_t xres, uint32_t yres, const T *data, bool doTri = false,
00047 float maxAniso = 8.f, ImageWrap wrapMode = TEXTURE_REPEAT);
00048 ~MIPMap();
00049 uint32_t Width() const { return width; }
00050 uint32_t Height() const { return height; }
00051 uint32_t Levels() const { return nLevels; }
00052 const T &Texel(uint32_t level, int s, int t) const;
00053 T Lookup(float s, float t, float width = 0.f) const;
00054 T Lookup(float s, float t, float ds0, float dt0,
00055 float ds1, float dt1) const;
00056 private:
00057
00058 struct ResampleWeight;
00059 ResampleWeight *resampleWeights(uint32_t oldres, uint32_t newres) {
00060 Assert(newres >= oldres);
00061 ResampleWeight *wt = new ResampleWeight[newres];
00062 float filterwidth = 2.f;
00063 for (uint32_t i = 0; i < newres; ++i) {
00064
00065 float center = (i + .5f) * oldres / newres;
00066 wt[i].firstTexel = Floor2Int((center - filterwidth) + 0.5f);
00067 for (int j = 0; j < 4; ++j) {
00068 float pos = wt[i].firstTexel + j + .5f;
00069 wt[i].weight[j] = Lanczos((pos - center) / filterwidth);
00070 }
00071
00072
00073 float invSumWts = 1.f / (wt[i].weight[0] + wt[i].weight[1] +
00074 wt[i].weight[2] + wt[i].weight[3]);
00075 for (uint32_t j = 0; j < 4; ++j)
00076 wt[i].weight[j] *= invSumWts;
00077 }
00078 return wt;
00079 }
00080 float clamp(float v) { return Clamp(v, 0.f, INFINITY); }
00081 RGBSpectrum clamp(const RGBSpectrum &v) { return v.Clamp(0.f, INFINITY); }
00082 SampledSpectrum clamp(const SampledSpectrum &v) { return v.Clamp(0.f, INFINITY); }
00083 T triangle(uint32_t level, float s, float t) const;
00084 T EWA(uint32_t level, float s, float t, float ds0, float dt0, float ds1, float dt1) const;
00085
00086
00087 bool doTrilinear;
00088 float maxAnisotropy;
00089 ImageWrap wrapMode;
00090 struct ResampleWeight {
00091 int firstTexel;
00092 float weight[4];
00093 };
00094 BlockedArray<T> **pyramid;
00095 uint32_t width, height, nLevels;
00096 #define WEIGHT_LUT_SIZE 128
00097 static float *weightLut;
00098 };
00099
00100
00101
00102
00103 template <typename T>
00104 MIPMap<T>::MIPMap(uint32_t sres, uint32_t tres, const T *img, bool doTri,
00105 float maxAniso, ImageWrap wm) {
00106 doTrilinear = doTri;
00107 maxAnisotropy = maxAniso;
00108 wrapMode = wm;
00109 T *resampledImage = NULL;
00110 if (!IsPowerOf2(sres) || !IsPowerOf2(tres)) {
00111
00112 uint32_t sPow2 = RoundUpPow2(sres), tPow2 = RoundUpPow2(tres);
00113
00114
00115 ResampleWeight *sWeights = resampleWeights(sres, sPow2);
00116 resampledImage = new T[sPow2 * tPow2];
00117
00118
00119 for (uint32_t t = 0; t < tres; ++t) {
00120 for (uint32_t s = 0; s < sPow2; ++s) {
00121
00122 resampledImage[t*sPow2+s] = 0.;
00123 for (int j = 0; j < 4; ++j) {
00124 int origS = sWeights[s].firstTexel + j;
00125 if (wrapMode == TEXTURE_REPEAT)
00126 origS = Mod(origS, sres);
00127 else if (wrapMode == TEXTURE_CLAMP)
00128 origS = Clamp(origS, 0, sres-1);
00129 if (origS >= 0 && origS < (int)sres)
00130 resampledImage[t*sPow2+s] += sWeights[s].weight[j] *
00131 img[t*sres + origS];
00132 }
00133 }
00134 }
00135 delete[] sWeights;
00136
00137
00138 ResampleWeight *tWeights = resampleWeights(tres, tPow2);
00139 T *workData = new T[tPow2];
00140 for (uint32_t s = 0; s < sPow2; ++s) {
00141 for (uint32_t t = 0; t < tPow2; ++t) {
00142 workData[t] = 0.;
00143 for (uint32_t j = 0; j < 4; ++j) {
00144 int offset = tWeights[t].firstTexel + j;
00145 if (wrapMode == TEXTURE_REPEAT) offset = Mod(offset, tres);
00146 else if (wrapMode == TEXTURE_CLAMP) offset = Clamp(offset, 0, tres-1);
00147 if (offset >= 0 && offset < (int)tres)
00148 workData[t] += tWeights[t].weight[j] *
00149 resampledImage[offset*sPow2 + s];
00150 }
00151 }
00152 for (uint32_t t = 0; t < tPow2; ++t)
00153 resampledImage[t*sPow2 + s] = clamp(workData[t]);
00154 }
00155 delete[] workData;
00156 delete[] tWeights;
00157 img = resampledImage;
00158 sres = sPow2;
00159 tres = tPow2;
00160 }
00161 width = sres;
00162 height = tres;
00163
00164 nLevels = 1 + Log2Int(float(max(sres, tres)));
00165 pyramid = new BlockedArray<T> *[nLevels];
00166
00167
00168 pyramid[0] = new BlockedArray<T>(sres, tres, img);
00169 for (uint32_t i = 1; i < nLevels; ++i) {
00170
00171 uint32_t sRes = max(1u, pyramid[i-1]->uSize()/2);
00172 uint32_t tRes = max(1u, pyramid[i-1]->vSize()/2);
00173 pyramid[i] = new BlockedArray<T>(sRes, tRes);
00174
00175
00176 for (uint32_t t = 0; t < tRes; ++t)
00177 for (uint32_t s = 0; s < sRes; ++s)
00178 (*pyramid[i])(s, t) = .25f *
00179 (Texel(i-1, 2*s, 2*t) + Texel(i-1, 2*s+1, 2*t) +
00180 Texel(i-1, 2*s, 2*t+1) + Texel(i-1, 2*s+1, 2*t+1));
00181 }
00182 if (resampledImage) delete[] resampledImage;
00183
00184 if (!weightLut) {
00185 weightLut = AllocAligned<float>(WEIGHT_LUT_SIZE);
00186 for (int i = 0; i < WEIGHT_LUT_SIZE; ++i) {
00187 float alpha = 2;
00188 float r2 = float(i) / float(WEIGHT_LUT_SIZE - 1);
00189 weightLut[i] = expf(-alpha * r2) - expf(-alpha);
00190 }
00191 }
00192 }
00193
00194
00195 template <typename T>
00196 const T &MIPMap<T>::Texel(uint32_t level, int s, int t) const {
00197 Assert(level < nLevels);
00198 const BlockedArray<T> &l = *pyramid[level];
00199
00200 switch (wrapMode) {
00201 case TEXTURE_REPEAT:
00202 s = Mod(s, l.uSize());
00203 t = Mod(t, l.vSize());
00204 break;
00205 case TEXTURE_CLAMP:
00206 s = Clamp(s, 0, l.uSize() - 1);
00207 t = Clamp(t, 0, l.vSize() - 1);
00208 break;
00209 case TEXTURE_BLACK: {
00210 static const T black = 0.f;
00211 if (s < 0 || s >= (int)l.uSize() ||
00212 t < 0 || t >= (int)l.vSize())
00213 return black;
00214 break;
00215 }
00216 }
00217 PBRT_ACCESSED_TEXEL(const_cast<MIPMap<T> *>(this), level, s, t);
00218 return l(s, t);
00219 }
00220
00221
00222 template <typename T>
00223 MIPMap<T>::~MIPMap() {
00224 for (uint32_t i = 0; i < nLevels; ++i)
00225 delete pyramid[i];
00226 delete[] pyramid;
00227 }
00228
00229
00230 template <typename T>
00231 T MIPMap<T>::Lookup(float s, float t, float width) const {
00232
00233 float level = nLevels - 1 + Log2(max(width, 1e-8f));
00234
00235
00236 PBRT_MIPMAP_TRILINEAR_FILTER(const_cast<MIPMap<T> *>(this), s, t, width, level, nLevels);
00237 if (level < 0)
00238 return triangle(0, s, t);
00239 else if (level >= nLevels - 1)
00240 return Texel(nLevels-1, 0, 0);
00241 else {
00242 uint32_t iLevel = Floor2Int(level);
00243 float delta = level - iLevel;
00244 return (1.f-delta) * triangle(iLevel, s, t) +
00245 delta * triangle(iLevel+1, s, t);
00246 }
00247 }
00248
00249
00250 template <typename T>
00251 T MIPMap<T>::triangle(uint32_t level, float s, float t) const {
00252 level = Clamp(level, 0, nLevels-1);
00253 s = s * pyramid[level]->uSize() - 0.5f;
00254 t = t * pyramid[level]->vSize() - 0.5f;
00255 int s0 = Floor2Int(s), t0 = Floor2Int(t);
00256 float ds = s - s0, dt = t - t0;
00257 return (1.f-ds) * (1.f-dt) * Texel(level, s0, t0) +
00258 (1.f-ds) * dt * Texel(level, s0, t0+1) +
00259 ds * (1.f-dt) * Texel(level, s0+1, t0) +
00260 ds * dt * Texel(level, s0+1, t0+1);
00261 }
00262
00263
00264 template <typename T>
00265 T MIPMap<T>::Lookup(float s, float t, float ds0, float dt0,
00266 float ds1, float dt1) const {
00267 if (doTrilinear) {
00268 PBRT_STARTED_TRILINEAR_TEXTURE_LOOKUP(s, t);
00269 T val = Lookup(s, t,
00270 2.f * max(max(fabsf(ds0), fabsf(dt0)),
00271 max(fabsf(ds1), fabsf(dt1))));
00272 PBRT_FINISHED_TRILINEAR_TEXTURE_LOOKUP();
00273 return val;
00274 }
00275 PBRT_STARTED_EWA_TEXTURE_LOOKUP(s, t);
00276
00277 if (ds0*ds0 + dt0*dt0 < ds1*ds1 + dt1*dt1) {
00278 swap(ds0, ds1);
00279 swap(dt0, dt1);
00280 }
00281 float majorLength = sqrtf(ds0*ds0 + dt0*dt0);
00282 float minorLength = sqrtf(ds1*ds1 + dt1*dt1);
00283
00284
00285 if (minorLength * maxAnisotropy < majorLength && minorLength > 0.f) {
00286 float scale = majorLength / (minorLength * maxAnisotropy);
00287 ds1 *= scale;
00288 dt1 *= scale;
00289 minorLength *= scale;
00290 }
00291 if (minorLength == 0.f) {
00292 PBRT_FINISHED_EWA_TEXTURE_LOOKUP();
00293 PBRT_STARTED_TRILINEAR_TEXTURE_LOOKUP(s, t);
00294 T val = triangle(0, s, t);
00295 PBRT_FINISHED_TRILINEAR_TEXTURE_LOOKUP();
00296 return val;
00297 }
00298
00299
00300 float lod = max(0.f, nLevels - 1.f + Log2(minorLength));
00301 uint32_t ilod = Floor2Int(lod);
00302 PBRT_MIPMAP_EWA_FILTER(const_cast<MIPMap<T> *>(this), s, t, ds0, ds1, dt0, dt1, minorLength, majorLength, lod, nLevels);
00303 float d = lod - ilod;
00304 T val = (1.f - d) * EWA(ilod, s, t, ds0, dt0, ds1, dt1) +
00305 d * EWA(ilod+1, s, t, ds0, dt0, ds1, dt1);
00306 PBRT_FINISHED_EWA_TEXTURE_LOOKUP();
00307 return val;
00308 }
00309
00310
00311 template <typename T>
00312 T MIPMap<T>::EWA(uint32_t level, float s, float t, float ds0, float dt0,
00313 float ds1, float dt1) const {
00314 if (level >= nLevels) return Texel(nLevels-1, 0, 0);
00315
00316 s = s * pyramid[level]->uSize() - 0.5f;
00317 t = t * pyramid[level]->vSize() - 0.5f;
00318 ds0 *= pyramid[level]->uSize();
00319 dt0 *= pyramid[level]->vSize();
00320 ds1 *= pyramid[level]->uSize();
00321 dt1 *= pyramid[level]->vSize();
00322
00323
00324 float A = dt0*dt0 + dt1*dt1 + 1;
00325 float B = -2.f * (ds0*dt0 + ds1*dt1);
00326 float C = ds0*ds0 + ds1*ds1 + 1;
00327 float invF = 1.f / (A*C - B*B*0.25f);
00328 A *= invF;
00329 B *= invF;
00330 C *= invF;
00331
00332
00333 float det = -B*B + 4.f*A*C;
00334 float invDet = 1.f / det;
00335 float uSqrt = sqrtf(det * C), vSqrt = sqrtf(A * det);
00336 int s0 = Ceil2Int (s - 2.f * invDet * uSqrt);
00337 int s1 = Floor2Int(s + 2.f * invDet * uSqrt);
00338 int t0 = Ceil2Int (t - 2.f * invDet * vSqrt);
00339 int t1 = Floor2Int(t + 2.f * invDet * vSqrt);
00340
00341
00342 T sum(0.);
00343 float sumWts = 0.f;
00344 for (int it = t0; it <= t1; ++it) {
00345 float tt = it - t;
00346 for (int is = s0; is <= s1; ++is) {
00347 float ss = is - s;
00348
00349 float r2 = A*ss*ss + B*ss*tt + C*tt*tt;
00350 if (r2 < 1.) {
00351 float weight = weightLut[min(Float2Int(r2 * WEIGHT_LUT_SIZE),
00352 WEIGHT_LUT_SIZE-1)];
00353 sum += Texel(level, is, it) * weight;
00354 sumWts += weight;
00355 }
00356 }
00357 }
00358 return sum / sumWts;
00359 }
00360
00361
00362 template <typename T> float *MIPMap<T>::weightLut = NULL;
00363
00364 #endif // PBRT_CORE_MIPMAP_H