|
|
// Copyright (c) 2009-2011 Ignacio Castano <castano@gmail.com>
|
|
|
//
|
|
|
// Permission is hereby granted, free of charge, to any person
|
|
|
// obtaining a copy of this software and associated documentation
|
|
|
// files (the "Software"), to deal in the Software without
|
|
|
// restriction, including without limitation the rights to use,
|
|
|
// copy, modify, merge, publish, distribute, sublicense, and/or sell
|
|
|
// copies of the Software, and to permit persons to whom the
|
|
|
// Software is furnished to do so, subject to the following
|
|
|
// conditions:
|
|
|
//
|
|
|
// The above copyright notice and this permission notice shall be
|
|
|
// included in all copies or substantial portions of the Software.
|
|
|
//
|
|
|
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
|
|
|
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
|
|
|
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
|
|
|
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
|
|
|
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
|
|
|
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
|
|
|
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
|
|
|
// OTHER DEALINGS IN THE SOFTWARE.
|
|
|
|
|
|
#include "CubeSurface.h"
|
|
|
#include "Surface.h"
|
|
|
|
|
|
#include "nvimage/DirectDrawSurface.h"
|
|
|
|
|
|
#include "nvmath/Vector.inl"
|
|
|
|
|
|
#include "nvcore/Array.inl"
|
|
|
#include "nvcore/StrLib.h"
|
|
|
|
|
|
using namespace nv;
|
|
|
using namespace nvtt;
|
|
|
|
|
|
|
|
|
|
|
|
// Solid angle of an axis aligned quad from (0,0,1) to (x,y,1)
|
|
|
// See: http://www.fizzmoll11.com/thesis/ for a derivation of this formula.
|
|
|
static float areaElement(float x, float y) {
|
|
|
return atan2f(x*y, sqrtf(x*x + y*y + 1));
|
|
|
}
|
|
|
|
|
|
// Solid angle of a hemicube texel.
|
|
|
static float solidAngleTerm(uint x, uint y, float inverseEdgeLength) {
|
|
|
// Transform x,y to [-1, 1] range, offset by 0.5 to point to texel center.
|
|
|
float u = (float(x) + 0.5f) * (2 * inverseEdgeLength) - 1.0f;
|
|
|
float v = (float(y) + 0.5f) * (2 * inverseEdgeLength) - 1.0f;
|
|
|
nvDebugCheck(u >= -1.0f && u <= 1.0f);
|
|
|
nvDebugCheck(v >= -1.0f && v <= 1.0f);
|
|
|
|
|
|
#if 1
|
|
|
// Exact solid angle:
|
|
|
float x0 = u - inverseEdgeLength;
|
|
|
float y0 = v - inverseEdgeLength;
|
|
|
float x1 = u + inverseEdgeLength;
|
|
|
float y1 = v + inverseEdgeLength;
|
|
|
float solidAngle = areaElement(x0, y0) - areaElement(x0, y1) - areaElement(x1, y0) + areaElement(x1, y1);
|
|
|
nvDebugCheck(solidAngle > 0.0f);
|
|
|
|
|
|
return solidAngle;
|
|
|
#else
|
|
|
// This formula is equivalent, but not as precise.
|
|
|
float pixel_area = nv::square(2.0f * inverseEdgeLength);
|
|
|
float dist_square = 1.0f + nv::square(u) + nv::square(v);
|
|
|
float cos_theta = 1.0f / sqrt(dist_square);
|
|
|
float cos_theta_d2 = cos_theta / dist_square; // Funny this is just 1/dist^3 or cos(tetha)^3
|
|
|
|
|
|
return pixel_area * cos_theta_d2;
|
|
|
#endif
|
|
|
}
|
|
|
|
|
|
|
|
|
static Vector3 texelDirection(uint face, uint x, uint y, int edgeLength, EdgeFixup fixupMethod)
|
|
|
{
|
|
|
float u, v;
|
|
|
if (fixupMethod == EdgeFixup_Stretch) {
|
|
|
// Transform x,y to [-1, 1] range, match up edges exactly.
|
|
|
u = float(x) * 2.0f / (edgeLength - 1) - 1.0f;
|
|
|
v = float(y) * 2.0f / (edgeLength - 1) - 1.0f;
|
|
|
}
|
|
|
else {
|
|
|
// Transform x,y to [-1, 1] range, offset by 0.5 to point to texel center.
|
|
|
u = (float(x) + 0.5f) * (2.0f / edgeLength) - 1.0f;
|
|
|
v = (float(y) + 0.5f) * (2.0f / edgeLength) - 1.0f;
|
|
|
}
|
|
|
|
|
|
if (fixupMethod == EdgeFixup_Warp) {
|
|
|
// Warp texel centers in the proximity of the edges.
|
|
|
float a = powf(float(edgeLength), 2.0f) / powf(float(edgeLength - 1), 3.0f);
|
|
|
u = a * powf(u, 3) + u;
|
|
|
v = a * powf(v, 3) + v;
|
|
|
}
|
|
|
|
|
|
nvDebugCheck(u >= -1.0f && u <= 1.0f);
|
|
|
nvDebugCheck(v >= -1.0f && v <= 1.0f);
|
|
|
|
|
|
Vector3 n;
|
|
|
|
|
|
if (face == 0) {
|
|
|
n.x = 1;
|
|
|
n.y = -v;
|
|
|
n.z = -u;
|
|
|
}
|
|
|
if (face == 1) {
|
|
|
n.x = -1;
|
|
|
n.y = -v;
|
|
|
n.z = u;
|
|
|
}
|
|
|
|
|
|
if (face == 2) {
|
|
|
n.x = u;
|
|
|
n.y = 1;
|
|
|
n.z = v;
|
|
|
}
|
|
|
if (face == 3) {
|
|
|
n.x = u;
|
|
|
n.y = -1;
|
|
|
n.z = -v;
|
|
|
}
|
|
|
|
|
|
if (face == 4) {
|
|
|
n.x = u;
|
|
|
n.y = -v;
|
|
|
n.z = 1;
|
|
|
}
|
|
|
if (face == 5) {
|
|
|
n.x = -u;
|
|
|
n.y = -v;
|
|
|
n.z = -1;
|
|
|
}
|
|
|
|
|
|
return normalizeFast(n);
|
|
|
}
|
|
|
|
|
|
|
|
|
TexelTable::TexelTable(uint edgeLength) : size(edgeLength) {
|
|
|
|
|
|
uint hsize = size/2;
|
|
|
|
|
|
// Allocate a small solid angle table that takes into account cube map symmetry.
|
|
|
solidAngleArray.resize(hsize * hsize);
|
|
|
|
|
|
for (uint y = 0; y < hsize; y++) {
|
|
|
for (uint x = 0; x < hsize; x++) {
|
|
|
solidAngleArray[y * hsize + x] = solidAngleTerm(hsize+x, hsize+y, 1.0f/edgeLength);
|
|
|
}
|
|
|
}
|
|
|
|
|
|
|
|
|
directionArray.resize(size*size*6);
|
|
|
|
|
|
for (uint f = 0; f < 6; f++) {
|
|
|
for (uint y = 0; y < size; y++) {
|
|
|
for (uint x = 0; x < size; x++) {
|
|
|
directionArray[(f * size + y) * size + x] = texelDirection(f, x, y, edgeLength, EdgeFixup_None);
|
|
|
}
|
|
|
}
|
|
|
}
|
|
|
}
|
|
|
|
|
|
const Vector3 & TexelTable::direction(uint f, uint x, uint y) const {
|
|
|
nvDebugCheck(f < 6 && x < size && y < size);
|
|
|
return directionArray[(f * size + y) * size + x];
|
|
|
}
|
|
|
|
|
|
float TexelTable::solidAngle(uint f, uint x, uint y) const {
|
|
|
uint hsize = size/2;
|
|
|
if (x >= hsize) x -= hsize;
|
|
|
else if (x < hsize) x = hsize - x - 1;
|
|
|
if (y >= hsize) y -= hsize;
|
|
|
else if (y < hsize) y = hsize - y - 1;
|
|
|
|
|
|
return solidAngleArray[y * hsize + x];
|
|
|
}
|
|
|
|
|
|
|
|
|
static const Vector3 faceNormals[6] = {
|
|
|
Vector3(1, 0, 0),
|
|
|
Vector3(-1, 0, 0),
|
|
|
Vector3(0, 1, 0),
|
|
|
Vector3(0, -1, 0),
|
|
|
Vector3(0, 0, 1),
|
|
|
Vector3(0, 0, -1),
|
|
|
};
|
|
|
|
|
|
static const Vector3 faceU[6] = {
|
|
|
Vector3(0, 0, -1),
|
|
|
Vector3(0, 0, 1),
|
|
|
Vector3(1, 0, 0),
|
|
|
Vector3(1, 0, 0),
|
|
|
Vector3(1, 0, 0),
|
|
|
Vector3(-1, 0, 0),
|
|
|
};
|
|
|
|
|
|
static const Vector3 faceV[6] = {
|
|
|
Vector3(0, -1, 0),
|
|
|
Vector3(0, -1, 0),
|
|
|
Vector3(0, 0, 1),
|
|
|
Vector3(0, 0, -1),
|
|
|
Vector3(0, -1, 0),
|
|
|
Vector3(0, -1, 0),
|
|
|
};
|
|
|
|
|
|
|
|
|
static Vector2 toPolar(Vector3::Arg v) {
|
|
|
Vector2 p;
|
|
|
p.x = atan2f(v.x, v.y); // theta
|
|
|
p.y = acosf(v.z); // phi
|
|
|
return p;
|
|
|
}
|
|
|
|
|
|
static Vector2 toPlane(float theta, float phi) {
|
|
|
float x = sinf(phi) * cosf(theta);
|
|
|
float y = sinf(phi) * sinf(theta);
|
|
|
float z = cosf(phi);
|
|
|
|
|
|
Vector2 p;
|
|
|
p.x = x / fabsf(z);
|
|
|
p.y = y / fabsf(z);
|
|
|
//p.x = tan(phi) * cos(theta);
|
|
|
//p.y = tan(phi) * sin(theta);
|
|
|
|
|
|
return p;
|
|
|
}
|
|
|
|
|
|
static Vector2 toPlane(Vector3::Arg v) {
|
|
|
Vector2 p;
|
|
|
p.x = v.x / fabsf(v.z);
|
|
|
p.y = v.y / fabsf(v.z);
|
|
|
return p;
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
CubeSurface::CubeSurface() : m(new CubeSurface::Private())
|
|
|
{
|
|
|
m->addRef();
|
|
|
}
|
|
|
|
|
|
CubeSurface::CubeSurface(const CubeSurface & cube) : m(cube.m)
|
|
|
{
|
|
|
if (m != NULL) m->addRef();
|
|
|
}
|
|
|
|
|
|
CubeSurface::~CubeSurface()
|
|
|
{
|
|
|
if (m != NULL) m->release();
|
|
|
m = NULL;
|
|
|
}
|
|
|
|
|
|
void CubeSurface::operator=(const CubeSurface & cube)
|
|
|
{
|
|
|
if (cube.m != NULL) cube.m->addRef();
|
|
|
if (m != NULL) m->release();
|
|
|
m = cube.m;
|
|
|
}
|
|
|
|
|
|
void CubeSurface::detach()
|
|
|
{
|
|
|
if (m->refCount() > 1)
|
|
|
{
|
|
|
m->release();
|
|
|
m = new CubeSurface::Private(*m);
|
|
|
m->addRef();
|
|
|
nvDebugCheck(m->refCount() == 1);
|
|
|
}
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
bool CubeSurface::isNull() const
|
|
|
{
|
|
|
return m->edgeLength == 0;
|
|
|
}
|
|
|
|
|
|
int CubeSurface::edgeLength() const
|
|
|
{
|
|
|
return m->edgeLength;
|
|
|
}
|
|
|
|
|
|
int CubeSurface::countMipmaps() const
|
|
|
{
|
|
|
return nv::countMipmaps(m->edgeLength);
|
|
|
}
|
|
|
|
|
|
Surface & CubeSurface::face(int f)
|
|
|
{
|
|
|
nvDebugCheck(f >= 0 && f < 6);
|
|
|
return m->face[f];
|
|
|
}
|
|
|
|
|
|
const Surface & CubeSurface::face(int f) const
|
|
|
{
|
|
|
nvDebugCheck(f >= 0 && f < 6);
|
|
|
return m->face[f];
|
|
|
}
|
|
|
|
|
|
|
|
|
bool CubeSurface::load(const char * fileName, int mipmap)
|
|
|
{
|
|
|
if (strEqual(Path::extension(fileName), ".dds")) {
|
|
|
nv::DirectDrawSurface dds;
|
|
|
if (!dds.load(fileName)) {
|
|
|
return false;
|
|
|
}
|
|
|
|
|
|
if (!dds.isValid()/* || !dds.isSupported()*/) {
|
|
|
return false;
|
|
|
}
|
|
|
|
|
|
if (!dds.isTextureCube()) {
|
|
|
return false;
|
|
|
}
|
|
|
|
|
|
// Make sure it's a valid cube.
|
|
|
if (dds.header.width != dds.header.height) return false;
|
|
|
//if ((dds.header.caps.caps2 & DDSCAPS2_CUBEMAP_ALL_FACES) != DDSCAPS2_CUBEMAP_ALL_FACES) return false;
|
|
|
|
|
|
if (mipmap < 0) {
|
|
|
mipmap = dds.mipmapCount() - 1 - mipmap;
|
|
|
}
|
|
|
if (mipmap < 0 || mipmap > I32(dds.mipmapCount())) return false;
|
|
|
|
|
|
|
|
|
nvtt::InputFormat inputFormat = nvtt::InputFormat_RGBA_16F;
|
|
|
|
|
|
if (dds.header.hasDX10Header()) {
|
|
|
if (dds.header.header10.dxgiFormat == DXGI_FORMAT_R16G16B16A16_FLOAT) inputFormat = nvtt::InputFormat_RGBA_16F;
|
|
|
else if (dds.header.header10.dxgiFormat == DXGI_FORMAT_R32G32B32A32_FLOAT) inputFormat = nvtt::InputFormat_RGBA_32F;
|
|
|
else if (dds.header.header10.dxgiFormat == DXGI_FORMAT_R32_FLOAT) inputFormat = nvtt::InputFormat_R_32F;
|
|
|
else return false;
|
|
|
}
|
|
|
else {
|
|
|
if ((dds.header.pf.flags & DDPF_FOURCC) != 0) {
|
|
|
if (dds.header.pf.fourcc == D3DFMT_A16B16G16R16F) inputFormat = nvtt::InputFormat_RGBA_16F;
|
|
|
else if (dds.header.pf.fourcc == D3DFMT_A32B32G32R32F) inputFormat = nvtt::InputFormat_RGBA_32F;
|
|
|
else if (dds.header.pf.fourcc == D3DFMT_R32F) inputFormat = nvtt::InputFormat_R_32F;
|
|
|
else return false;
|
|
|
}
|
|
|
else {
|
|
|
if (dds.header.pf.bitcount == 32 /*&& ...*/) inputFormat = nvtt::InputFormat_BGRA_8UB;
|
|
|
else return false; // @@ Do pixel format conversions!
|
|
|
}
|
|
|
}
|
|
|
|
|
|
uint edgeLength = dds.surfaceWidth(mipmap);
|
|
|
uint size = dds.surfaceSize(mipmap);
|
|
|
|
|
|
void * data = malloc(size);
|
|
|
|
|
|
for (int f = 0; f < 6; f++) {
|
|
|
dds.readSurface(f, mipmap, data, size);
|
|
|
m->face[f].setImage(inputFormat, edgeLength, edgeLength, 1, data);
|
|
|
}
|
|
|
|
|
|
m->edgeLength = edgeLength;
|
|
|
|
|
|
free(data);
|
|
|
|
|
|
return true;
|
|
|
}
|
|
|
|
|
|
return false;
|
|
|
}
|
|
|
|
|
|
bool CubeSurface::save(const char * fileName) const
|
|
|
{
|
|
|
// @@ TODO
|
|
|
return false;
|
|
|
}
|
|
|
|
|
|
struct ivec2 {
|
|
|
uint x;
|
|
|
uint y;
|
|
|
};
|
|
|
// posx negx posy negy posz negz
|
|
|
static const ivec2 foldOffsetVerticalCross[6] = { {2, 1}, {0, 1}, {1, 0}, {1, 2}, {1, 1}, {1, 3} };
|
|
|
static const ivec2 foldOffsetHorizontalCross[6] = { {2, 1}, {0, 1}, {1, 0}, {1, 2}, {1, 1}, {3, 1} };
|
|
|
static const ivec2 foldOffsetColumn[6] = { {0, 0}, {0, 1}, {0, 2}, {0, 3}, {0, 4}, {0, 5} };
|
|
|
static const ivec2 foldOffsetRow[6] = { {0, 0}, {1, 0}, {2, 0}, {3, 0}, {4, 0}, {5, 0} };
|
|
|
|
|
|
void CubeSurface::fold(const Surface & tex, CubeLayout layout)
|
|
|
{
|
|
|
ivec2 const* offsets = 0;
|
|
|
uint edgeLength;
|
|
|
|
|
|
switch(layout) {
|
|
|
case CubeLayout_LatitudeLongitude:
|
|
|
case CubeLayout_VerticalCross:
|
|
|
edgeLength = tex.height() / 4;
|
|
|
offsets = foldOffsetVerticalCross;
|
|
|
break;
|
|
|
case CubeLayout_HorizontalCross:
|
|
|
edgeLength = tex.width() / 4;
|
|
|
offsets = foldOffsetHorizontalCross;
|
|
|
break;
|
|
|
case CubeLayout_Column:
|
|
|
edgeLength = tex.width();
|
|
|
offsets = foldOffsetColumn;
|
|
|
break;
|
|
|
case CubeLayout_Row:
|
|
|
edgeLength = tex.height();
|
|
|
offsets = foldOffsetRow;
|
|
|
break;
|
|
|
}
|
|
|
|
|
|
m->edgeLength = edgeLength;
|
|
|
for(uint f = 0; f < 6; f++) {
|
|
|
uint x = offsets[f].x * edgeLength;
|
|
|
uint y = offsets[f].y * edgeLength;
|
|
|
m->face[f] = tex.createSubImage(x, x + edgeLength - 1, y, y + edgeLength - 1, 0, 0);
|
|
|
}
|
|
|
|
|
|
if(layout == CubeLayout_VerticalCross || layout == CubeLayout_LatitudeLongitude) {
|
|
|
// Back face needs to be rotated 180 degrees
|
|
|
m->face[5].flipX();
|
|
|
m->face[5].flipY();
|
|
|
}
|
|
|
}
|
|
|
|
|
|
Surface CubeSurface::unfold(CubeLayout layout) const
|
|
|
{
|
|
|
ivec2 const* offsets = 0;
|
|
|
uint edgeLength = m->edgeLength;
|
|
|
uint width;
|
|
|
uint height;
|
|
|
|
|
|
switch(layout) {
|
|
|
case CubeLayout_LatitudeLongitude:
|
|
|
case CubeLayout_VerticalCross:
|
|
|
offsets = foldOffsetVerticalCross;
|
|
|
width = 3 * edgeLength;
|
|
|
height = 4 * edgeLength;
|
|
|
// Back face needs to be rotated 180 degrees
|
|
|
m->face[5].flipX();
|
|
|
m->face[5].flipY();
|
|
|
break;
|
|
|
case CubeLayout_HorizontalCross:
|
|
|
offsets = foldOffsetHorizontalCross;
|
|
|
width = 4 * edgeLength;
|
|
|
height = 3 * edgeLength;
|
|
|
break;
|
|
|
case CubeLayout_Column:
|
|
|
offsets = foldOffsetColumn;
|
|
|
width = edgeLength;
|
|
|
height = 6 * edgeLength;
|
|
|
break;
|
|
|
case CubeLayout_Row:
|
|
|
offsets = foldOffsetRow;
|
|
|
width = 6 * edgeLength;
|
|
|
height = edgeLength;
|
|
|
break;
|
|
|
}
|
|
|
|
|
|
Surface surface;
|
|
|
surface.setImage(width, height, 1);
|
|
|
for(uint f = 0; f < 6; f++) {
|
|
|
uint x = offsets[f].x * edgeLength;
|
|
|
uint y = offsets[f].y * edgeLength;
|
|
|
surface.copy(m->face[f], 0, 0, 0, edgeLength, edgeLength, 1, x, y, 0);
|
|
|
}
|
|
|
|
|
|
if(layout == CubeLayout_VerticalCross || layout == CubeLayout_LatitudeLongitude) {
|
|
|
// Undo back face rotation
|
|
|
m->face[5].flipY();
|
|
|
m->face[5].flipX();
|
|
|
}
|
|
|
return surface;
|
|
|
}
|
|
|
|
|
|
float CubeSurface::average(int channel) const
|
|
|
{
|
|
|
const uint edgeLength = m->edgeLength;
|
|
|
m->allocateTexelTable();
|
|
|
|
|
|
float total = 0.0f;
|
|
|
float sum = 0.0f;
|
|
|
|
|
|
for (int f = 0; f < 6; f++) {
|
|
|
float * c = m->face[f].m->image->channel(channel);
|
|
|
|
|
|
for (uint y = 0; y < edgeLength; y++) {
|
|
|
for (uint x = 0; x < edgeLength; x++) {
|
|
|
float solidAngle = m->texelTable->solidAngle(f, x, y);
|
|
|
|
|
|
total += solidAngle;
|
|
|
sum += c[y * edgeLength + x] * solidAngle;
|
|
|
}
|
|
|
}
|
|
|
}
|
|
|
|
|
|
return sum / total;
|
|
|
}
|
|
|
|
|
|
void CubeSurface::range(int channel, float * minimum_ptr, float * maximum_ptr) const
|
|
|
{
|
|
|
const uint edgeLength = m->edgeLength;
|
|
|
m->allocateTexelTable();
|
|
|
|
|
|
float minimum = NV_FLOAT_MAX;
|
|
|
float maximum = 0.0f;
|
|
|
|
|
|
for (int f = 0; f < 6; f++) {
|
|
|
float * c = m->face[f].m->image->channel(channel);
|
|
|
|
|
|
for (uint y = 0; y < edgeLength; y++) {
|
|
|
for (uint x = 0; x < edgeLength; x++) {
|
|
|
|
|
|
minimum = nv::min(minimum, c[y * edgeLength + x]);
|
|
|
maximum = nv::max(maximum, c[y * edgeLength + x]);
|
|
|
}
|
|
|
}
|
|
|
}
|
|
|
|
|
|
*minimum_ptr = minimum;
|
|
|
*maximum_ptr = maximum;
|
|
|
}
|
|
|
|
|
|
void CubeSurface::clamp(int channel, float low/*= 0.0f*/, float high/*= 1.0f*/) {
|
|
|
for (int f = 0; f < 6; f++) {
|
|
|
m->face[f].clamp(channel, low, high);
|
|
|
}
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
#include "nvmath/SphericalHarmonic.h"
|
|
|
|
|
|
CubeSurface CubeSurface::irradianceFilter(int size, EdgeFixup fixupMethod) const
|
|
|
{
|
|
|
// Evaluate spherical harmonic for each output texel.
|
|
|
CubeSurface output;
|
|
|
output.m->allocate(size);
|
|
|
|
|
|
Sh2 shr, shg, shb;
|
|
|
|
|
|
computeIrradianceSH3(0, shr.coef);
|
|
|
computeIrradianceSH3(1, shg.coef);
|
|
|
computeIrradianceSH3(2, shb.coef);
|
|
|
|
|
|
// @@ Sample spherical harmonic from every direction.
|
|
|
|
|
|
return CubeSurface();
|
|
|
}
|
|
|
|
|
|
|
|
|
void CubeSurface::computeLuminanceIrradianceSH3(float coef[9]) const{
|
|
|
|
|
|
m->allocateTexelTable();
|
|
|
|
|
|
// Transform this cube to spherical harmonic basis
|
|
|
Sh2 sh;
|
|
|
|
|
|
// For each texel of the input cube.
|
|
|
const uint edgeLength = m->edgeLength;
|
|
|
for (uint f = 0; f < 6; f++) {
|
|
|
|
|
|
const Surface & inputFace = m->face[f];
|
|
|
const FloatImage * inputImage = inputFace.m->image;
|
|
|
|
|
|
for (uint y = 0; y < edgeLength; y++) {
|
|
|
for (uint x = 0; x < edgeLength; x++) {
|
|
|
|
|
|
Vector3 dir = m->texelTable->direction(f, x, y);
|
|
|
float solidAngle = m->texelTable->solidAngle(f, x, y);
|
|
|
|
|
|
Sh2 shDir;
|
|
|
shDir.eval(dir);
|
|
|
|
|
|
float r = inputImage->pixel(0, x, y, 0);
|
|
|
float g = inputImage->pixel(1, x, y, 0);
|
|
|
float b = inputImage->pixel(2, x, y, 0);
|
|
|
float lum = 0.333f * (r + g + b); // @@ use the proper luminance formula.
|
|
|
|
|
|
sh.addScaled(shDir, lum * solidAngle);
|
|
|
}
|
|
|
}
|
|
|
}
|
|
|
|
|
|
for (int i = 0; i < 9; i++) {
|
|
|
coef[i] = sh.coef[i];
|
|
|
}
|
|
|
}
|
|
|
|
|
|
|
|
|
void CubeSurface::computeIrradianceSH3(int channel, float coef[9]) const {
|
|
|
|
|
|
m->allocateTexelTable();
|
|
|
|
|
|
// Transform this cube to spherical harmonic basis
|
|
|
Sh2 sh;
|
|
|
|
|
|
// For each texel of the input cube.
|
|
|
const uint edgeLength = m->edgeLength;
|
|
|
for (uint f = 0; f < 6; f++) {
|
|
|
|
|
|
const Surface & inputFace = m->face[f];
|
|
|
const FloatImage * inputImage = inputFace.m->image;
|
|
|
|
|
|
for (uint y = 0; y < edgeLength; y++) {
|
|
|
for (uint x = 0; x < edgeLength; x++) {
|
|
|
|
|
|
Vector3 dir = m->texelTable->direction(f, x, y);
|
|
|
float solidAngle = m->texelTable->solidAngle(f, x, y);
|
|
|
|
|
|
Sh2 shDir;
|
|
|
shDir.eval(dir);
|
|
|
|
|
|
float c = inputImage->pixel(channel, x, y, 0);
|
|
|
|
|
|
sh.addScaled(shDir, c * solidAngle);
|
|
|
}
|
|
|
}
|
|
|
}
|
|
|
|
|
|
for (int i = 0; i < 9; i++) {
|
|
|
coef[i] = sh.elemAt(i);
|
|
|
}
|
|
|
}
|
|
|
|
|
|
|
|
|
// Convolve filter against this cube.
|
|
|
Vector3 CubeSurface::Private::applyAngularFilter(const Vector3 & filterDir, float coneAngle, float * filterTable, int tableSize)
|
|
|
{
|
|
|
const float cosineConeAngle = cosf(coneAngle);
|
|
|
nvDebugCheck(cosineConeAngle >= 0);
|
|
|
|
|
|
Vector3 color(0);
|
|
|
float sum = 0;
|
|
|
|
|
|
// Things I have tried to speed this up:
|
|
|
// - Compute accurate bounds assuming cone axis aligned to plane, result was too small elsewhere.
|
|
|
// - Compute ellipse that results in the cone/plane intersection and compute its bounds. Sometimes intersection is a parabolla, hard to handle that case.
|
|
|
// - Compute the 6 axis aligned planes that bound the cone, clip faces against planes. Resulting plane equations are way too complex.
|
|
|
|
|
|
// What AMD CubeMapGen does:
|
|
|
// - Compute conservative bounds on the primary face, wrap around the adjacent faces.
|
|
|
|
|
|
|
|
|
// For each texel of the input cube.
|
|
|
for (uint f = 0; f < 6; f++) {
|
|
|
|
|
|
// Test face cone agains filter cone.
|
|
|
float cosineFaceAngle = dot(filterDir, faceNormals[f]);
|
|
|
float faceAngle = acosf(cosineFaceAngle);
|
|
|
|
|
|
if (faceAngle > coneAngle + atanf(sqrtf(2))) {
|
|
|
// Skip face.
|
|
|
continue;
|
|
|
}
|
|
|
|
|
|
const int L = I32(edgeLength-1);
|
|
|
int x0 = 0, x1 = L;
|
|
|
int y0 = 0, y1 = L;
|
|
|
|
|
|
#if 0
|
|
|
float u0 = -1;
|
|
|
float u1 = 1;
|
|
|
float v0 = -1;
|
|
|
float v1 = 1;
|
|
|
|
|
|
// @@ Compute uvs.
|
|
|
|
|
|
// Expand uv coordinates from [-1,1] to [0, edgeLength)
|
|
|
u0 = (u0 + 1) * edgeLength * 0.5f - 0.5f;
|
|
|
v0 = (v0 + 1) * edgeLength * 0.5f - 0.5f;
|
|
|
u1 = (u1 + 1) * edgeLength * 0.5f - 0.5f;
|
|
|
v1 = (v1 + 1) * edgeLength * 0.5f - 0.5f;
|
|
|
nvDebugCheck(u0 >= -0.5f && u0 <= edgeLength - 0.5f);
|
|
|
nvDebugCheck(v0 >= -0.5f && v0 <= edgeLength - 0.5f);
|
|
|
nvDebugCheck(u1 >= -0.5f && u1 <= edgeLength - 0.5f);
|
|
|
nvDebugCheck(v1 >= -0.5f && v1 <= edgeLength - 0.5f);
|
|
|
|
|
|
x0 = clamp(ifloor(u0), 0, L);
|
|
|
y0 = clamp(ifloor(v0), 0, L);
|
|
|
x1 = clamp(iceil(u1), 0, L);
|
|
|
y1 = clamp(iceil(v1), 0, L);
|
|
|
#endif
|
|
|
|
|
|
nvDebugCheck(x1 >= x0);
|
|
|
nvDebugCheck(y1 >= y0);
|
|
|
|
|
|
if (x1 == x0 || y1 == y0) {
|
|
|
// Skip this face.
|
|
|
continue;
|
|
|
}
|
|
|
|
|
|
|
|
|
const Surface & inputFace = face[f];
|
|
|
const FloatImage * inputImage = inputFace.m->image;
|
|
|
|
|
|
for (int y = y0; y <= y1; y++) {
|
|
|
bool inside = false;
|
|
|
for (int x = x0; x <= x1; x++) {
|
|
|
|
|
|
Vector3 dir = texelTable->direction(f, x, y);
|
|
|
float cosineAngle = dot(dir, filterDir);
|
|
|
|
|
|
if (cosineAngle > cosineConeAngle) {
|
|
|
float solidAngle = texelTable->solidAngle(f, x, y);
|
|
|
//float scale = powf(saturate(cosineAngle), cosinePower);
|
|
|
|
|
|
int idx = int(saturate(cosineAngle) * (tableSize - 1));
|
|
|
float scale = filterTable[idx]; // @@ Do bilinear interpolation?
|
|
|
|
|
|
float contribution = solidAngle * scale;
|
|
|
|
|
|
sum += contribution;
|
|
|
color.x += contribution * inputImage->pixel(0, x, y, 0);
|
|
|
color.y += contribution * inputImage->pixel(1, x, y, 0);
|
|
|
color.z += contribution * inputImage->pixel(2, x, y, 0);
|
|
|
|
|
|
inside = true;
|
|
|
}
|
|
|
else if (inside) {
|
|
|
// Filter scale is monotonic, if we have been inside once and we just exit, then we can skip the rest of the row.
|
|
|
// We could do the same thing for the columns and skip entire rows.
|
|
|
break;
|
|
|
}
|
|
|
}
|
|
|
}
|
|
|
}
|
|
|
|
|
|
color *= (1.0f / sum);
|
|
|
|
|
|
return color;
|
|
|
}
|
|
|
|
|
|
// We want to find the alpha such that:
|
|
|
// cos(alpha)^cosinePower = epsilon
|
|
|
// That's: acos(epsilon^(1/cosinePower))
|
|
|
|
|
|
// We can cull texels in two different ways:
|
|
|
// - culling faces that do not touch the cone.
|
|
|
// - computing one rectangle per face, find intersection between cone and face.
|
|
|
// -
|
|
|
|
|
|
// Other speedups:
|
|
|
// - parallelize. Done.
|
|
|
// - use ISPC?
|
|
|
|
|
|
|
|
|
// Convolve filter against this cube.
|
|
|
Vector3 CubeSurface::Private::applyCosinePowerFilter(const Vector3 & filterDir, float coneAngle, float cosinePower)
|
|
|
{
|
|
|
const float cosineConeAngle = cosf(coneAngle);
|
|
|
nvDebugCheck(cosineConeAngle >= 0);
|
|
|
|
|
|
Vector3 color(0);
|
|
|
float sum = 0;
|
|
|
|
|
|
// Things I have tried to speed this up:
|
|
|
// - Compute accurate bounds assuming cone axis aligned to plane, result was too small elsewhere.
|
|
|
// - Compute ellipse that results in the cone/plane intersection and compute its bounds. Sometimes intersection is a parabolla, hard to handle that case.
|
|
|
// - Compute the 6 axis aligned planes that bound the cone, clip faces against planes. Resulting plane equations are way too complex.
|
|
|
|
|
|
// What AMD CubeMapGen does:
|
|
|
// - Compute conservative bounds on the primary face, wrap around the adjacent faces.
|
|
|
|
|
|
|
|
|
// For each texel of the input cube.
|
|
|
for (uint f = 0; f < 6; f++) {
|
|
|
|
|
|
// Test face cone agains filter cone.
|
|
|
float cosineFaceAngle = dot(filterDir, faceNormals[f]);
|
|
|
float faceAngle = acosf(cosineFaceAngle);
|
|
|
|
|
|
if (faceAngle > coneAngle + atanf(sqrtf(2))) {
|
|
|
// Skip face.
|
|
|
continue;
|
|
|
}
|
|
|
|
|
|
const int L = I32(edgeLength-1);
|
|
|
int x0 = 0, x1 = L;
|
|
|
int y0 = 0, y1 = L;
|
|
|
|
|
|
#if 0
|
|
|
float u0 = -1;
|
|
|
float u1 = 1;
|
|
|
float v0 = -1;
|
|
|
float v1 = 1;
|
|
|
|
|
|
// @@ Compute uvs.
|
|
|
|
|
|
// Expand uv coordinates from [-1,1] to [0, edgeLength)
|
|
|
u0 = (u0 + 1) * edgeLength * 0.5f - 0.5f;
|
|
|
v0 = (v0 + 1) * edgeLength * 0.5f - 0.5f;
|
|
|
u1 = (u1 + 1) * edgeLength * 0.5f - 0.5f;
|
|
|
v1 = (v1 + 1) * edgeLength * 0.5f - 0.5f;
|
|
|
nvDebugCheck(u0 >= -0.5f && u0 <= edgeLength - 0.5f);
|
|
|
nvDebugCheck(v0 >= -0.5f && v0 <= edgeLength - 0.5f);
|
|
|
nvDebugCheck(u1 >= -0.5f && u1 <= edgeLength - 0.5f);
|
|
|
nvDebugCheck(v1 >= -0.5f && v1 <= edgeLength - 0.5f);
|
|
|
|
|
|
x0 = clamp(ifloor(u0), 0, L);
|
|
|
y0 = clamp(ifloor(v0), 0, L);
|
|
|
x1 = clamp(iceil(u1), 0, L);
|
|
|
y1 = clamp(iceil(v1), 0, L);
|
|
|
#endif
|
|
|
|
|
|
nvDebugCheck(x1 >= x0);
|
|
|
nvDebugCheck(y1 >= y0);
|
|
|
|
|
|
if (x1 == x0 || y1 == y0) {
|
|
|
// Skip this face.
|
|
|
continue;
|
|
|
}
|
|
|
|
|
|
|
|
|
const Surface & inputFace = face[f];
|
|
|
const FloatImage * inputImage = inputFace.m->image;
|
|
|
|
|
|
for (int y = y0; y <= y1; y++) {
|
|
|
bool inside = false;
|
|
|
for (int x = x0; x <= x1; x++) {
|
|
|
|
|
|
Vector3 dir = texelTable->direction(f, x, y);
|
|
|
float cosineAngle = dot(dir, filterDir);
|
|
|
|
|
|
if (cosineAngle > cosineConeAngle) {
|
|
|
float solidAngle = texelTable->solidAngle(f, x, y);
|
|
|
float scale = powf(saturate(cosineAngle), cosinePower);
|
|
|
float contribution = solidAngle * scale;
|
|
|
|
|
|
sum += contribution;
|
|
|
color.x += contribution * inputImage->pixel(0, x, y, 0);
|
|
|
color.y += contribution * inputImage->pixel(1, x, y, 0);
|
|
|
color.z += contribution * inputImage->pixel(2, x, y, 0);
|
|
|
|
|
|
inside = true;
|
|
|
}
|
|
|
else if (inside) {
|
|
|
// Filter scale is monotonic, if we have been inside once and we just exit, then we can skip the rest of the row.
|
|
|
// We could do the same thing for the columns and skip entire rows.
|
|
|
break;
|
|
|
}
|
|
|
}
|
|
|
}
|
|
|
}
|
|
|
|
|
|
color *= (1.0f / sum);
|
|
|
|
|
|
return color;
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
#include "nvthread/ParallelFor.h"
|
|
|
|
|
|
struct ApplyAngularFilterContext {
|
|
|
CubeSurface::Private * inputCube;
|
|
|
CubeSurface::Private * filteredCube;
|
|
|
float coneAngle;
|
|
|
float * filterTable;
|
|
|
int tableSize;
|
|
|
EdgeFixup fixupMethod;
|
|
|
};
|
|
|
|
|
|
void ApplyAngularFilterTask(void * context, int id)
|
|
|
{
|
|
|
ApplyAngularFilterContext * ctx = (ApplyAngularFilterContext *)context;
|
|
|
|
|
|
int size = ctx->filteredCube->edgeLength;
|
|
|
|
|
|
int f = id / (size * size);
|
|
|
int idx = id % (size * size);
|
|
|
int y = idx / size;
|
|
|
int x = idx % size;
|
|
|
|
|
|
nvtt::Surface & filteredFace = ctx->filteredCube->face[f];
|
|
|
FloatImage * filteredImage = filteredFace.m->image;
|
|
|
|
|
|
const Vector3 filterDir = texelDirection(f, x, y, size, ctx->fixupMethod);
|
|
|
|
|
|
// Convolve filter against cube.
|
|
|
Vector3 color = ctx->inputCube->applyAngularFilter(filterDir, ctx->coneAngle, ctx->filterTable, ctx->tableSize);
|
|
|
|
|
|
filteredImage->pixel(0, idx) = color.x;
|
|
|
filteredImage->pixel(1, idx) = color.y;
|
|
|
filteredImage->pixel(2, idx) = color.z;
|
|
|
}
|
|
|
|
|
|
|
|
|
CubeSurface CubeSurface::cosinePowerFilter(int size, float cosinePower, EdgeFixup fixupMethod) const
|
|
|
{
|
|
|
// Allocate output cube.
|
|
|
CubeSurface filteredCube;
|
|
|
filteredCube.m->allocate(size);
|
|
|
|
|
|
// Texel table is stored along with the surface so that it's computed only once.
|
|
|
m->allocateTexelTable();
|
|
|
|
|
|
const float threshold = 0.001f;
|
|
|
const float coneAngle = acosf(powf(threshold, 1.0f/cosinePower));
|
|
|
|
|
|
|
|
|
// For each texel of the output cube.
|
|
|
/*for (uint f = 0; f < 6; f++) {
|
|
|
nvtt::Surface filteredFace = filteredCube.m->face[f];
|
|
|
FloatImage * filteredImage = filteredFace.m->image;
|
|
|
|
|
|
for (uint y = 0; y < uint(size); y++) {
|
|
|
for (uint x = 0; x < uint(size); x++) {
|
|
|
|
|
|
const Vector3 filterDir = texelDirection(f, x, y, size, fixupMethod);
|
|
|
|
|
|
// Convolve filter against cube.
|
|
|
Vector3 color = m->applyCosinePowerFilter(filterDir, coneAngle, cosinePower);
|
|
|
|
|
|
filteredImage->pixel(0, x, y, 0) = color.x;
|
|
|
filteredImage->pixel(1, x, y, 0) = color.y;
|
|
|
filteredImage->pixel(2, x, y, 0) = color.z;
|
|
|
}
|
|
|
}
|
|
|
}*/
|
|
|
|
|
|
ApplyAngularFilterContext context;
|
|
|
context.inputCube = m;
|
|
|
context.filteredCube = filteredCube.m;
|
|
|
context.coneAngle = coneAngle;
|
|
|
context.fixupMethod = fixupMethod;
|
|
|
|
|
|
context.tableSize = 512;
|
|
|
context.filterTable = new float[context.tableSize];
|
|
|
|
|
|
// @@ Instead of looking up table between [0 - 1] we should probably use [cos(coneAngle), 1]
|
|
|
|
|
|
for (int i = 0; i < context.tableSize; i++) {
|
|
|
float f = float(i) / (context.tableSize - 1);
|
|
|
context.filterTable[i] = powf(f, cosinePower);
|
|
|
}
|
|
|
|
|
|
|
|
|
nv::ParallelFor parallelFor(ApplyAngularFilterTask, &context);
|
|
|
parallelFor.run(6 * size * size);
|
|
|
|
|
|
// @@ Implement edge averaging.
|
|
|
if (fixupMethod == EdgeFixup_Average) {
|
|
|
for (uint f = 0; f < 6; f++) {
|
|
|
nvtt::Surface filteredFace = filteredCube.m->face[f];
|
|
|
FloatImage * filteredImage = filteredFace.m->image;
|
|
|
|
|
|
// For each component.
|
|
|
for (uint c = 0; c < 3; c++) {
|
|
|
// @@ For each corner, sample the two adjacent faces.
|
|
|
filteredImage->pixel(c, 0, 0, 0);
|
|
|
filteredImage->pixel(c, size-1, 0, 0);
|
|
|
filteredImage->pixel(c, 0, size-1, 0);
|
|
|
filteredImage->pixel(c, size-1, size-1, 0);
|
|
|
|
|
|
// @@ For each edge, sample the adjacent face.
|
|
|
|
|
|
}
|
|
|
}
|
|
|
}
|
|
|
|
|
|
return filteredCube;
|
|
|
}
|
|
|
|
|
|
|
|
|
// Sample cubemap in the given direction.
|
|
|
Vector3 CubeSurface::Private::sample(const Vector3 & dir)
|
|
|
{
|
|
|
int f = -1;
|
|
|
if (fabs(dir.x) > fabs(dir.y) && fabs(dir.x) > fabs(dir.z)) {
|
|
|
if (dir.x > 0) f = 0;
|
|
|
else f = 1;
|
|
|
}
|
|
|
else if (fabs(dir.y) > fabs(dir.z)) {
|
|
|
if (dir.y > 0) f = 2;
|
|
|
else f = 3;
|
|
|
}
|
|
|
else {
|
|
|
if (dir.z > 0) f = 4;
|
|
|
else f = 5;
|
|
|
}
|
|
|
nvDebugCheck(f != -1);
|
|
|
|
|
|
// uv coordinates corresponding to filterDir.
|
|
|
float u = dot(dir, faceU[f]);
|
|
|
float v = dot(dir, faceV[f]);
|
|
|
|
|
|
FloatImage * img = face[f].m->image;
|
|
|
|
|
|
Vector3 color;
|
|
|
color.x = img->sampleLinearClamp(0, u, v);
|
|
|
color.y = img->sampleLinearClamp(1, u, v);
|
|
|
color.z = img->sampleLinearClamp(2, u, v);
|
|
|
|
|
|
return color;
|
|
|
}
|
|
|
|
|
|
// @@ Not tested!
|
|
|
CubeSurface CubeSurface::fastResample(int size, EdgeFixup fixupMethod) const
|
|
|
{
|
|
|
// Allocate output cube.
|
|
|
CubeSurface resampledCube;
|
|
|
resampledCube.m->allocate(size);
|
|
|
|
|
|
// For each texel of the output cube.
|
|
|
for (uint f = 0; f < 6; f++) {
|
|
|
nvtt::Surface resampledFace = resampledCube.m->face[f];
|
|
|
FloatImage * resampledImage = resampledFace.m->image;
|
|
|
|
|
|
for (uint y = 0; y < uint(size); y++) {
|
|
|
for (uint x = 0; x < uint(size); x++) {
|
|
|
|
|
|
const Vector3 filterDir = texelDirection(f, x, y, size, fixupMethod);
|
|
|
|
|
|
Vector3 color = m->sample(filterDir);
|
|
|
|
|
|
resampledImage->pixel(0, x, y, 0) = color.x;
|
|
|
resampledImage->pixel(1, x, y, 0) = color.y;
|
|
|
resampledImage->pixel(2, x, y, 0) = color.z;
|
|
|
}
|
|
|
}
|
|
|
}
|
|
|
|
|
|
// @@ Implement edge averaging. Share this code with cosinePowerFilter
|
|
|
if (fixupMethod == EdgeFixup_Average) {
|
|
|
}
|
|
|
|
|
|
return resampledCube;
|
|
|
}
|
|
|
|
|
|
void CubeSurface::toLinear(float gamma)
|
|
|
{
|
|
|
if (isNull()) return;
|
|
|
|
|
|
detach();
|
|
|
|
|
|
for (int i = 0; i < 6; i++) {
|
|
|
m->face[i].toLinear(gamma);
|
|
|
}
|
|
|
}
|
|
|
|
|
|
void CubeSurface::toGamma(float gamma)
|
|
|
{
|
|
|
if (isNull()) return;
|
|
|
|
|
|
detach();
|
|
|
|
|
|
for (int i = 0; i < 6; i++) {
|
|
|
m->face[i].toGamma(gamma);
|
|
|
}
|
|
|
}
|
|
|
|
|
|
|
|
|
#if 0
|
|
|
// @@ Provide solar azimuth.
|
|
|
#include "ArHoseSkyModel.h"
|
|
|
void CubeSurface::sky(float turbidity, float albedo[3], float solarElevation) {
|
|
|
|
|
|
ArHosekSkyModelState * skymodel_state[3];
|
|
|
|
|
|
for (int i = 0; i < num_channels; i++) {
|
|
|
skymodel_state[i] = arhosekskymodelstate_alloc_init(turbidity, albedo[i], solarElevation);
|
|
|
}
|
|
|
|
|
|
// 700 nm (red), 546.1 nm (green) and 435.8 nm (blue).
|
|
|
float channel_center[3] = {
|
|
|
700, // Red 620<32>740,
|
|
|
546.1, // Green 520<32>570,
|
|
|
435.8, // Blue 450<35>490,
|
|
|
};
|
|
|
|
|
|
// @@ For each pixel:
|
|
|
// What's the channel center for the RGB model?
|
|
|
double skydome_result[3];
|
|
|
for (unsigned int i = 0; i < num_channels; i++) {
|
|
|
skydome_result[i] = arhosekskymodel_radiance(skymodel_state[i], theta, gamma, channel_center[i]);
|
|
|
}
|
|
|
|
|
|
for (int i = 0; i < num_channels; i++) {
|
|
|
arhosek_skymodelstate_free(skymodel_state[i]);
|
|
|
}
|
|
|
|
|
|
/*
|
|
|
ArHosekXYZSkyModelState * skymodel_state[3];
|
|
|
|
|
|
for (int i = 0; i < num_channels; i++) {
|
|
|
skymodel_state[i] = arhosek_xyz_skymodelstate_alloc_init(turbidity, albedo[i], solarElevation);
|
|
|
}
|
|
|
|
|
|
// @@ For each pixel.
|
|
|
double skydome_result[3];
|
|
|
for (unsigned int i = 0; i < num_channels; i++) {
|
|
|
skydome_result[i] = arhosek_xyz_skymodel_radiance(skymodel_state[i], theta, gamma, i);
|
|
|
}
|
|
|
|
|
|
for (int i = 0; i < num_channels; i++) {
|
|
|
arhosek_xyz_skymodelstate_free(skymodel_state[i]);
|
|
|
}
|
|
|
*/
|
|
|
}
|
|
|
#endif
|