Go to the documentation of this file.
18 #ifndef OPENVDB_TOOLS_RAYTRACER_HAS_BEEN_INCLUDED
19 #define OPENVDB_TOOLS_RAYTRACER_HAS_BEEN_INCLUDED
33 #include <type_traits>
36 #ifdef OPENVDB_TOOLS_RAYTRACER_USE_EXR
37 #include <OpenEXR/ImfPixelType.h>
38 #include <OpenEXR/ImfChannelList.h>
39 #include <OpenEXR/ImfOutputFile.h>
40 #include <OpenEXR/ImfHeader.h>
41 #include <OpenEXR/ImfFrameBuffer.h>
54 template<
typename Gr
idT>
58 size_t pixelSamples = 1,
59 unsigned int seed = 0,
60 bool threaded =
true);
63 template<
typename Gr
idT,
typename IntersectorT>
68 size_t pixelSamples = 1,
69 unsigned int seed = 0,
70 bool threaded =
true);
77 template<
typename Gr
idT,
typename IntersectorT = tools::LevelSetRayIntersector<Gr
idT> >
82 using Vec3Type =
typename IntersectorT::Vec3Type;
83 using RayType =
typename IntersectorT::RayType;
89 size_t pixelSamples = 1,
90 unsigned int seed = 0);
97 size_t pixelSamples = 1,
98 unsigned int seed = 0);
107 void setGrid(
const GridT& grid);
111 void setIntersector(
const IntersectorT& inter);
129 void setPixelSamples(
size_t pixelSamples,
unsigned int seed = 0);
132 void render(
bool threaded =
true)
const;
136 void operator()(
const tbb::blocked_range<size_t>& range)
const;
139 const bool mIsMaster;
142 std::unique_ptr<const BaseShader> mShader;
154 template <
typename IntersectorT,
typename SamplerT = tools::BoxSampler>
160 using RayType =
typename IntersectorT::RayType;
164 static_assert(std::is_floating_point<ValueType>::value,
165 "VolumeRender requires a floating-point-valued grid");
174 void render(
bool threaded=
true)
const;
181 void setIntersector(
const IntersectorT& inter);
213 void print(std::ostream& os = std::cout,
int verboseLevel = 1);
217 void operator()(
const tbb::blocked_range<size_t>& range)
const;
221 AccessorType mAccessor;
223 std::unique_ptr<IntersectorT> mPrimary, mShadow;
224 Real mPrimaryStep, mShadowStep, mCutOff, mLightGain;
225 Vec3R mLightDir, mLightColor, mAbsorption, mScattering;
241 RGBA() : r(0), g(0), b(0), a(1) {}
242 explicit RGBA(
ValueT intensity) : r(intensity), g(intensity), b(intensity), a(1) {}
244 r(_r), g(_g), b(_b), a(_a)
246 RGBA(
double _r,
double _g,
double _b,
double _a = 1.0)
247 : r(static_cast<
ValueT>(_r))
248 , g(static_cast<
ValueT>(_g))
249 , b(static_cast<
ValueT>(_b))
250 , a(static_cast<
ValueT>(_a))
260 const float s = rhs.
a*(1.0f-a);
271 Film(
size_t width,
size_t height)
272 : mWidth(width), mHeight(height), mSize(width*height), mPixels(new
RGBA[mSize])
276 : mWidth(width), mHeight(height), mSize(width*height), mPixels(new
RGBA[mSize])
285 return mPixels[w + h*mWidth];
292 return mPixels[w + h*mWidth];
295 void fill(
const RGBA& rgb=
RGBA(0)) {
for (
size_t i=0; i<mSize; ++i) mPixels[i] = rgb; }
298 RGBA *p = mPixels.get();
299 for (
size_t j = 0; j < mHeight; ++j) {
300 for (
size_t i = 0; i < mWidth; ++i, ++p) {
301 *p = ((i & size) ^ (j & size)) ? c1 : c2;
308 std::string
name(fileName);
309 if (
name.find_last_of(
".") == std::string::npos)
name.append(
".ppm");
311 std::unique_ptr<unsigned char[]> buffer(
new unsigned char[3*mSize]);
312 unsigned char *tmp = buffer.get(), *q = tmp;
313 RGBA* p = mPixels.get();
316 *q++ =
static_cast<unsigned char>(255.0f*(*p ).r);
317 *q++ =
static_cast<unsigned char>(255.0f*(*p ).g);
318 *q++ =
static_cast<unsigned char>(255.0f*(*p++).b);
321 std::ofstream os(
name.c_str(), std::ios_base::binary);
323 std::cerr <<
"Error opening PPM file \"" <<
name <<
"\"" << std::endl;
327 os <<
"P6\n" << mWidth <<
" " << mHeight <<
"\n255\n";
328 os.write(
reinterpret_cast<const char*
>(&(*tmp)), 3 * mSize *
sizeof(
unsigned char));
331 #ifdef OPENVDB_TOOLS_RAYTRACER_USE_EXR
332 void saveEXR(
const std::string& fileName,
size_t compression = 2,
size_t threads = 8)
334 std::string
name(fileName);
335 if (
name.find_last_of(
".") == std::string::npos)
name.append(
".exr");
337 if (threads>0) Imf::setGlobalThreadCount(threads);
338 Imf::Header header(mWidth, mHeight);
339 if (compression==0) header.compression() = Imf::NO_COMPRESSION;
340 if (compression==1) header.compression() = Imf::RLE_COMPRESSION;
341 if (compression>=2) header.compression() = Imf::ZIP_COMPRESSION;
342 header.channels().insert(
"R", Imf::Channel(Imf::FLOAT));
343 header.channels().insert(
"G", Imf::Channel(Imf::FLOAT));
344 header.channels().insert(
"B", Imf::Channel(Imf::FLOAT));
345 header.channels().insert(
"A", Imf::Channel(Imf::FLOAT));
347 Imf::FrameBuffer framebuffer;
348 framebuffer.insert(
"R", Imf::Slice( Imf::FLOAT, (
char *) &(mPixels[0].r),
349 sizeof (RGBA),
sizeof (RGBA) * mWidth));
350 framebuffer.insert(
"G", Imf::Slice( Imf::FLOAT, (
char *) &(mPixels[0].g),
351 sizeof (RGBA),
sizeof (RGBA) * mWidth));
352 framebuffer.insert(
"B", Imf::Slice( Imf::FLOAT, (
char *) &(mPixels[0].b),
353 sizeof (RGBA),
sizeof (RGBA) * mWidth));
354 framebuffer.insert(
"A", Imf::Slice( Imf::FLOAT, (
char *) &(mPixels[0].a),
355 sizeof (RGBA),
sizeof (RGBA) * mWidth));
357 Imf::OutputFile file(
name.c_str(), header);
358 file.setFrameBuffer(framebuffer);
359 file.writePixels(mHeight);
363 size_t width()
const {
return mWidth; }
364 size_t height()
const {
return mHeight; }
369 size_t mWidth, mHeight, mSize;
370 std::unique_ptr<RGBA[]> mPixels;
381 double frameWidth,
double nearPlane,
double farPlane)
383 , mScaleWidth(frameWidth)
384 , mScaleHeight(frameWidth * double(film.height()) / double(film.width()))
386 assert(nearPlane > 0 && farPlane > nearPlane);
390 mScreenToWorld.accumPostTranslation(translation);
391 this->initRay(nearPlane, farPlane);
398 size_t width()
const {
return mFilm->width(); }
399 size_t height()
const {
return mFilm->height(); }
407 const Vec3R orig = mScreenToWorld.applyMap(
Vec3R(0.0));
408 const Vec3R dir = orig - xyz;
410 Mat4d xform = math::aim<Mat4d>(dir, up);
411 xform.postTranslate(orig);
413 this->initRay(mRay.t0(), mRay.t1());
419 return Vec3R( (2 * i /
double(mFilm->width()) - 1) * mScaleWidth,
420 (1 - 2 * j /
double(mFilm->height())) * mScaleHeight, z );
427 size_t i,
size_t j,
double iOffset = 0.5,
double jOffset = 0.5)
const = 0;
432 mRay.setTimes(t0, t1);
433 mRay.setEye(mScreenToWorld.applyMap(
Vec3R(0.0)));
434 mRay.setDir(mScreenToWorld.applyJacobian(
Vec3R(0.0, 0.0, -1.0)));
465 double focalLength = 50.0,
466 double aperture = 41.2136,
467 double nearPlane = 1e-3,
469 :
BaseCamera(film,
rotation, translation, 0.5*aperture/focalLength, nearPlane, farPlane)
479 size_t i,
size_t j,
double iOffset = 0.5,
double jOffset = 0.5)
const override
482 Vec3R dir = BaseCamera::rasterToScreen(
Real(i) + iOffset,
Real(j) + jOffset, -1.0);
483 dir = BaseCamera::mScreenToWorld.applyJacobian(dir);
494 return 360.0 / M_PI * atan(aperture/(2.0*length));
500 return aperture/(2.0*(tan(fov * M_PI / 360.0)));
523 double frameWidth = 1.0,
524 double nearPlane = 1e-3,
532 size_t i,
size_t j,
double iOffset = 0.5,
double jOffset = 0.5)
const override
535 Vec3R eye = BaseCamera::rasterToScreen(
Real(i) + iOffset,
Real(j) + jOffset, 0.0);
536 ray.
setEye(BaseCamera::mScreenToWorld.applyMap(eye));
573 MatteShader(
const GridT& grid) : mAcc(grid.getAccessor()), mXform(&grid.transform()) {}
578 typename GridT::ValueType v = zeroVal<typename GridT::ValueType>();
579 SamplerType::sample(mAcc, mXform->worldToIndex(xyz), v);
585 typename GridT::ConstAccessor mAcc;
590 template<
typename SamplerType>
615 template<
typename GridT = Film::RGBA,
620 NormalShader(
const GridT& grid) : mAcc(grid.getAccessor()), mXform(&grid.transform()) {}
625 typename GridT::ValueType v = zeroVal<typename GridT::ValueType>();
626 SamplerType::sample(mAcc, mXform->worldToIndex(xyz), v);
627 return Film::RGBA(v[0]*(normal[0]+1.0), v[1]*(normal[1]+1.0), v[2]*(normal[2]+1.0));
632 typename GridT::ConstAccessor mAcc;
637 template<
typename SamplerType>
646 return mRGBA *
Film::RGBA(normal[0] + 1.0, normal[1] + 1.0, normal[2] + 1.0);
662 template<
typename GridT = Film::RGBA,
669 , mInvDim(1.0/bbox.extents())
670 , mAcc(grid.getAccessor())
671 , mXform(&grid.transform())
678 typename GridT::ValueType v = zeroVal<typename GridT::ValueType>();
679 SamplerType::sample(mAcc, mXform->worldToIndex(xyz), v);
680 const Vec3R rgb = (xyz - mMin) * mInvDim;
686 const Vec3R mMin, mInvDim;
687 typename GridT::ConstAccessor mAcc;
692 template<
typename SamplerType>
697 : mMin(bbox.
min()), mInvDim(1.0/bbox.extents()), mRGBA(c) {}
702 const Vec3R rgb = (xyz - mMin)*mInvDim;
703 return mRGBA*
Film::RGBA(rgb[0], rgb[1], rgb[2]);
708 const Vec3R mMin, mInvDim;
722 template<
typename GridT = Film::RGBA,
727 DiffuseShader(
const GridT& grid): mAcc(grid.getAccessor()), mXform(&grid.transform()) {}
732 typename GridT::ValueType v = zeroVal<typename GridT::ValueType>();
733 SamplerType::sample(mAcc, mXform->worldToIndex(xyz), v);
742 typename GridT::ConstAccessor mAcc;
747 template <
typename SamplerType>
776 template<
typename Gr
idT>
785 tracer(grid, shader, camera, pixelSamples, seed);
790 template<
typename Gr
idT,
typename IntersectorT>
792 const IntersectorT& inter,
807 template<
typename Gr
idT,
typename IntersectorT>
817 mShader(shader.copy()),
823 template<
typename Gr
idT,
typename IntersectorT>
833 mShader(shader.copy()),
839 template<
typename Gr
idT,
typename IntersectorT>
844 mInter(other.mInter),
845 mShader(other.mShader->copy()),
846 mCamera(other.mCamera),
847 mSubPixels(other.mSubPixels)
851 template<
typename Gr
idT,
typename IntersectorT>
855 if (mIsMaster)
delete [] mRand;
858 template<
typename Gr
idT,
typename IntersectorT>
863 mInter = IntersectorT(grid);
866 template<
typename Gr
idT,
typename IntersectorT>
874 template<
typename Gr
idT,
typename IntersectorT>
879 mShader.reset(shader.
copy());
882 template<
typename Gr
idT,
typename IntersectorT>
890 template<
typename Gr
idT,
typename IntersectorT>
895 if (pixelSamples == 0) {
898 mSubPixels = pixelSamples - 1;
900 if (mSubPixels > 0) {
901 mRand =
new double[16];
903 for (
size_t i=0; i<16; ++i) mRand[i] = rand();
909 template<
typename Gr
idT,
typename IntersectorT>
911 render(
bool threaded)
const
913 tbb::blocked_range<size_t> range(0, mCamera->height());
914 threaded ? tbb::parallel_for(range, *
this) : (*this)(range);
917 template<
typename Gr
idT,
typename IntersectorT>
919 operator()(
const tbb::blocked_range<size_t>& range)
const
923 const float frac = 1.0f / (1.0f + float(mSubPixels));
924 for (
size_t j=range.begin(), n=0, je = range.end(); j<je; ++j) {
925 for (
size_t i=0, ie = mCamera->width(); i<ie; ++i) {
927 RayType ray = mCamera->getRay(i, j);
928 Film::RGBA c = mInter.intersectsWS(ray, xyz, nml) ? shader(xyz, nml, ray.dir()) : bg;
929 for (
size_t k=0; k<mSubPixels; ++k, n +=2 ) {
930 ray = mCamera->getRay(i, j, mRand[n & 15], mRand[(n+1) & 15]);
931 c += mInter.intersectsWS(ray, xyz, nml) ? shader(xyz, nml, ray.dir()) : bg;
940 template<
typename IntersectorT,
typename SampleT>
943 : mAccessor(inter.grid().getConstAccessor())
945 , mPrimary(new IntersectorT(inter))
946 , mShadow(new IntersectorT(inter))
952 , mLightColor(0.7, 0.7, 0.7)
958 template<
typename IntersectorT,
typename SampleT>
961 : mAccessor(other.mAccessor)
962 , mCamera(other.mCamera)
963 , mPrimary(new IntersectorT(*(other.mPrimary)))
964 , mShadow(new IntersectorT(*(other.mShadow)))
965 , mPrimaryStep(other.mPrimaryStep)
966 , mShadowStep(other.mShadowStep)
967 , mCutOff(other.mCutOff)
968 , mLightGain(other.mLightGain)
969 , mLightDir(other.mLightDir)
970 , mLightColor(other.mLightColor)
971 , mAbsorption(other.mAbsorption)
972 , mScattering(other.mScattering)
976 template<
typename IntersectorT,
typename SampleT>
978 print(std::ostream& os,
int verboseLevel)
980 if (verboseLevel>0) {
981 os <<
"\nPrimary step: " << mPrimaryStep
982 <<
"\nShadow step: " << mShadowStep
983 <<
"\nCutoff: " << mCutOff
984 <<
"\nLightGain: " << mLightGain
985 <<
"\nLightDir: " << mLightDir
986 <<
"\nLightColor: " << mLightColor
987 <<
"\nAbsorption: " << mAbsorption
988 <<
"\nScattering: " << mScattering << std::endl;
990 mPrimary->print(os, verboseLevel);
993 template<
typename IntersectorT,
typename SampleT>
997 mPrimary.reset(
new IntersectorT(inter));
998 mShadow.reset(
new IntersectorT(inter));
1001 template<
typename IntersectorT,
typename SampleT>
1003 render(
bool threaded)
const
1005 tbb::blocked_range<size_t> range(0, mCamera->height());
1006 threaded ? tbb::parallel_for(range, *
this) : (*this)(range);
1009 template<
typename IntersectorT,
typename SampleT>
1011 operator()(
const tbb::blocked_range<size_t>& range)
const
1013 SamplerType sampler(mAccessor, mShadow->grid().transform());
1016 const Vec3R extinction = -mScattering-mAbsorption, One(1.0);
1017 const Vec3R albedo = mLightColor*mScattering/(mScattering+mAbsorption);
1018 const Real sGain = mLightGain;
1019 const Real pStep = mPrimaryStep;
1020 const Real sStep = mShadowStep;
1021 const Real cutoff = mCutOff;
1030 std::vector<typename RayType::TimeSpan> pTS, sTS;
1035 for (
size_t j=range.begin(), je = range.end(); j<je; ++j) {
1036 for (
size_t i=0, ie = mCamera->width(); i<ie; ++i) {
1038 bg.
a = bg.
r = bg.
g = bg.
b = 0;
1039 RayType pRay = mCamera->getRay(i, j);
1040 if( !mPrimary->setWorldRay(pRay))
continue;
1041 Vec3R pTrans(1.0), pLumi(0.0);
1044 while (mPrimary->march(pT0, pT1)) {
1045 for (
Real pT = pStep*ceil(pT0/pStep); pT <= pT1; pT += pStep) {
1047 mPrimary->hits(pTS);
1048 for (
size_t k=0; k<pTS.size(); ++k) {
1049 Real pT = pStep*ceil(pTS[k].t0/pStep), pT1=pTS[k].t1;
1050 for (; pT <= pT1; pT += pStep) {
1052 Vec3R pPos = mPrimary->getWorldPos(pT);
1054 if (density < cutoff)
continue;
1058 if( !mShadow->setWorldRay(sRay))
continue;
1061 while (mShadow->march(sT0, sT1)) {
1062 for (
Real sT = sStep*ceil(sT0/sStep); sT <= sT1; sT+= sStep) {
1065 for (
size_t l=0; l<sTS.size(); ++l) {
1066 Real sT = sStep*ceil(sTS[l].t0/sStep), sT1=sTS[l].t1;
1067 for (; sT <= sT1; sT+= sStep) {
1069 const Real d = sampler.
wsSample(mShadow->getWorldPos(sT));
1070 if (d < cutoff)
continue;
1071 sTrans *=
math::Exp(extinction * d * sStep/(1.0+sT*sGain));
1072 if (sTrans.
lengthSqr()<cutoff)
goto Luminance;
1076 pLumi += albedo * sTrans * pTrans * (One-dT);
1078 if (pTrans.lengthSqr()<cutoff)
goto Pixel;
1094 #endif // OPENVDB_TOOLS_RAYTRACER_HAS_BEEN_INCLUDED
@ Z_AXIS
Definition: Math.h:872
@ X_AXIS
Definition: Math.h:870
Simple generator of random numbers over the range [0, 1)
Definition: Math.h:129
MatType scale(const Vec3< typename MatType::value_type > &s)
Return a matrix that scales by s.
Definition: Mat.h:620
Coord Abs(const Coord &xyz)
Definition: Coord.h:515
@ Y_AXIS
Definition: Math.h:871
Accelerated intersection of a ray with a narrow-band level set or a generic (e.g. density) volume....
General-purpose arithmetic and comparison routines, most of which accept arbitrary value types (or at...
T lengthSqr() const
Definition: Vec3.h:209
Vec3< typename promote< T, typename Coord::ValueType >::type > operator+(const Vec3< T > &v0, const Coord &v1)
Allow a Coord to be added to or subtracted from a Vec3.
Definition: Coord.h:525
math::Vec3< Real > Vec3R
Definition: Types.h:49
void setDir(const Vec3Type &dir)
Definition: Ray.h:66
void setEye(const Vec3Type &eye)
Definition: Ray.h:64
Type Exp(const Type &x)
Return ex.
Definition: Math.h:677
MatType unit(const MatType &mat, typename MatType::value_type eps=1.0e-8)
Return a copy of the given matrix with its upper 3×3 rows normalized.
Definition: Mat.h:653
Vec3< T > unit(T eps=0) const
return normalized this, throws if null vector
Definition: Vec3.h:372
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:146
void scaleTimes(RealT scale)
Definition: Ray.h:84
Axis-aligned bounding box.
Definition: BBox.h:24
const Vec3T & dir() const
Definition: Ray.h:99
Mat4< double > Mat4d
Definition: Mat4.h:1334
const Name & name
Definition: PointAttribute.h:544
#define OPENVDB_VERSION_NAME
The version namespace name for this library version.
Definition: version.h:94
A general linear transform using homogeneous coordinates to perform rotation, scaling,...
Definition: Maps.h:299
T dot(const Vec3< T > &v) const
Dot product.
Definition: Vec3.h:189
bool normalize(T eps=T(1.0e-7))
this = normalized this
Definition: Vec3.h:360
MatType rotation(const Quat< typename MatType::value_type > &q, typename MatType::value_type eps=static_cast< typename MatType::value_type >(1.0e-8))
Return the rotation matrix specified by the given quaternion.
Definition: Mat.h:177
Definition: Exceptions.h:65
Definition: Exceptions.h:13
double Real
Definition: Types.h:37
Mat3< typename promote< T0, T1 >::type > operator*(const Mat3< T0 > &m0, const Mat3< T1 > &m1)
Multiply m0 by m1 and return the resulting matrix.
Definition: Mat3.h:618
#define OPENVDB_THROW(exception, message)
Definition: Exceptions.h:82