Improving the LowPassFilter to take advantage of the features made available by the voxel refactoring.

Also improved the test case.
This commit is contained in:
David Williams 2012-09-28 18:09:21 +02:00
parent bc01223237
commit 01e2a88b13
5 changed files with 145 additions and 104 deletions

View File

@ -101,12 +101,14 @@ int main(int argc, char *argv[])
createCubeInVolume(volData, Vector3DInt32(midPos-10, 1, midPos-10), Vector3DInt32(midPos+10, maxPos-1, midPos+10), MaterialDensityPair44::getMaxDensity());
createCubeInVolume(volData, Vector3DInt32(midPos-10, midPos-10 ,1), Vector3DInt32(midPos+10, midPos+10, maxPos-1), MaterialDensityPair44::getMaxDensity());
//Smooth part of the volume
RawVolume<MaterialDensityPair44> tempVolume(PolyVox::Region(0,0,0,128, 128, 128));
LowPassFilter< LargeVolume<MaterialDensityPair44>, RawVolume<MaterialDensityPair44> > pass1(&volData, PolyVox::Region(Vector3DInt32(62, 62, 62), Vector3DInt32(126, 126, 126)), &tempVolume, PolyVox::Region(Vector3DInt32(62, 62, 62), Vector3DInt32(126, 126, 126)), 3);
pass1.executeSAT();
LowPassFilter< RawVolume<MaterialDensityPair44>, LargeVolume<MaterialDensityPair44> > pass2(&tempVolume, PolyVox::Region(Vector3DInt32(62, 62, 62), Vector3DInt32(126, 126, 126)), &volData, PolyVox::Region(Vector3DInt32(62, 62, 62), Vector3DInt32(126, 126, 126)), 3);
pass2.executeSAT();
//I've removed this smoothing because it doesn't really make sense to apply a low pass filter to a volume with material values.
//I could implement the mathematical operators for MaterialDensityPair in such a way that they ignores the materials but this
//seems to be setting a bad example. Users can add this operators in their own classes if they want smoothing.
//RawVolume<MaterialDensityPair44> tempVolume(PolyVox::Region(0,0,0,128, 128, 128));
//LowPassFilter< LargeVolume<MaterialDensityPair44>, RawVolume<MaterialDensityPair44> > pass1(&volData, PolyVox::Region(Vector3DInt32(62, 62, 62), Vector3DInt32(126, 126, 126)), &tempVolume, PolyVox::Region(Vector3DInt32(62, 62, 62), Vector3DInt32(126, 126, 126)), 3);
//pass1.executeSAT();
//LowPassFilter< RawVolume<MaterialDensityPair44>, LargeVolume<MaterialDensityPair44> > pass2(&tempVolume, PolyVox::Region(Vector3DInt32(62, 62, 62), Vector3DInt32(126, 126, 126)), &volData, PolyVox::Region(Vector3DInt32(62, 62, 62), Vector3DInt32(126, 126, 126)), 3);
//pass2.executeSAT();
QApplication app(argc, argv);

View File

@ -36,17 +36,9 @@ freely, subject to the following restrictions:
namespace PolyVox
{
///This class represents a voxel storing only a density.
/// This class represents a voxel storing only a density.
////////////////////////////////////////////////////////////////////////////////
/// In order to perform a surface extraction on a LargeVolume, PolyVox needs the underlying
/// voxel type to provide both getDensity() and getMaterial() functions. The getDensity()
/// function is used to determine if a voxel is 'solid', and if it is then the getMaterial()
/// funtion is used to determine what material should be assigned to the resulting mesh.
///
/// This class meets these requirements, although it only actually stores a density value.
/// For the getMaterial() function it just returns a constant value of '1'.
///
/// \sa Material, MaterialDensityPair
/// Detailed description...
////////////////////////////////////////////////////////////////////////////////
// int32_t template parameter is a dummy, required as the compiler expects to be able to declare an
@ -62,9 +54,18 @@ namespace PolyVox
typedef Type DensityType;
typedef int32_t MaterialType; //Shouldn't define this one...
/// Constructor
Density() : m_uDensity(0) {}
/// Copy constructor
Density(DensityType uDensity) : m_uDensity(uDensity) {}
/// Copy constructor with cast
template <typename CastType> explicit Density(const Density<CastType>& density) throw()
{
m_uDensity = static_cast<Type>(density.getDensity());
}
bool operator==(const Density& rhs) const throw()
{
return (m_uDensity == rhs.m_uDensity);
@ -75,12 +76,22 @@ namespace PolyVox
return !(*this == rhs);
}
// For densities we can supply mathematical operators which behave in an intuitive way.
// In particular the ability to add and subtract densities is important in order to
// apply an averaging filter. The ability to divide by an integer is also needed for
// this same purpose.
Density<Type>& operator+=(const Density<Type>& rhs)
{
m_uDensity += rhs.m_uDensity;
return *this;
}
Density<Type>& operator-=(const Density<Type>& rhs)
{
m_uDensity -= rhs.m_uDensity;
return *this;
}
Density<Type>& operator/=(uint32_t rhs)
{
m_uDensity /= rhs;
@ -97,9 +108,36 @@ namespace PolyVox
DensityType m_uDensity;
};
template <typename Type>
Density<Type> operator+(const Density<Type>& lhs, const Density<Type>& rhs) throw()
{
Density<Type> result = lhs;
result += rhs;
return result;
}
template <typename Type>
Density<Type> operator-(const Density<Type>& lhs, const Density<Type>& rhs) throw()
{
Density<Type> result = lhs;
result -= rhs;
return result;
}
template <typename Type>
Density<Type> operator/(const Density<Type>& lhs, uint32_t rhs) throw()
{
Density<Type> result = lhs;
result /= rhs;
return result;
}
// These are the predefined density types. The 8-bit types are sufficient for many purposes (including
// most games) but 16-bit and float types do have uses particularly in medical/scientific visualisation.
typedef Density<uint8_t> Density8;
typedef Density<uint16_t> Density16;
typedef Density<uint32_t> Density32;
typedef Density<float> DensityFloat;
/**
* This is a specialisation of DefaultMarchingCubesController for the Density voxel type

View File

@ -30,7 +30,7 @@ freely, subject to the following restrictions:
namespace PolyVox
{
template< typename SrcVolumeType, typename DstVolumeType>
template< typename SrcVolumeType, typename DstVolumeType, typename AccumulationType>
class LowPassFilter
{
public:

View File

@ -23,8 +23,8 @@ freely, subject to the following restrictions:
namespace PolyVox
{
template< typename SrcVolumeType, typename DstVolumeType>
LowPassFilter<SrcVolumeType, DstVolumeType>::LowPassFilter(SrcVolumeType* pVolSrc, Region regSrc, DstVolumeType* pVolDst, Region regDst, uint32_t uKernelSize)
template< typename SrcVolumeType, typename DstVolumeType, typename AccumulationType>
LowPassFilter<SrcVolumeType, DstVolumeType, AccumulationType>::LowPassFilter(SrcVolumeType* pVolSrc, Region regSrc, DstVolumeType* pVolDst, Region regDst, uint32_t uKernelSize)
:m_pVolSrc(pVolSrc)
,m_regSrc(regSrc)
,m_pVolDst(pVolDst)
@ -43,8 +43,8 @@ namespace PolyVox
}
}
template< typename SrcVolumeType, typename DstVolumeType>
void LowPassFilter<SrcVolumeType, DstVolumeType>::execute()
template< typename SrcVolumeType, typename DstVolumeType, typename AccumulationType>
void LowPassFilter<SrcVolumeType, DstVolumeType, AccumulationType>::execute()
{
int32_t iSrcMinX = m_regSrc.getLowerCorner().getX();
int32_t iSrcMinY = m_regSrc.getLowerCorner().getY();
@ -70,52 +70,50 @@ namespace PolyVox
{
for(int32_t iSrcX = iSrcMinX, iDstX = iDstMinX; iSrcX <= iSrcMaxX; iSrcX++, iDstX++)
{
//VoxelType tSrcVoxel = m_pVolSrc->getVoxelAt(iSrcX, iSrcY, iSrcZ);
AccumulationType tSrcVoxel(0);
srcSampler.setPosition(iSrcX, iSrcY, iSrcZ);
typename SrcVolumeType::VoxelType tSrcVoxel = srcSampler.getVoxel();
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1nx1ny1nz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1nx1ny0pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1nx1ny1pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1nx0py1nz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1nx0py0pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1nx0py1pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1nx1py1nz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1nx1py0pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1nx1py1pz());
tSrcVoxel += srcSampler.peekVoxel1nx1ny1nz();
tSrcVoxel += srcSampler.peekVoxel1nx1ny0pz();
tSrcVoxel += srcSampler.peekVoxel1nx1ny1pz();
tSrcVoxel += srcSampler.peekVoxel1nx0py1nz();
tSrcVoxel += srcSampler.peekVoxel1nx0py0pz();
tSrcVoxel += srcSampler.peekVoxel1nx0py1pz();
tSrcVoxel += srcSampler.peekVoxel1nx1py1nz();
tSrcVoxel += srcSampler.peekVoxel1nx1py0pz();
tSrcVoxel += srcSampler.peekVoxel1nx1py1pz();
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel0px1ny1nz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel0px1ny0pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel0px1ny1pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel0px0py1nz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel0px0py0pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel0px0py1pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel0px1py1nz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel0px1py0pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel0px1py1pz());
tSrcVoxel += srcSampler.peekVoxel0px1ny1nz();
tSrcVoxel += srcSampler.peekVoxel0px1ny0pz();
tSrcVoxel += srcSampler.peekVoxel0px1ny1pz();
tSrcVoxel += srcSampler.peekVoxel0px0py1nz();
//tSrcVoxel += srcSampler.peekVoxel0px0py0pz();
tSrcVoxel += srcSampler.peekVoxel0px0py1pz();
tSrcVoxel += srcSampler.peekVoxel0px1py1nz();
tSrcVoxel += srcSampler.peekVoxel0px1py0pz();
tSrcVoxel += srcSampler.peekVoxel0px1py1pz();
tSrcVoxel += srcSampler.peekVoxel1px1ny1nz();
tSrcVoxel += srcSampler.peekVoxel1px1ny0pz();
tSrcVoxel += srcSampler.peekVoxel1px1ny1pz();
tSrcVoxel += srcSampler.peekVoxel1px0py1nz();
tSrcVoxel += srcSampler.peekVoxel1px0py0pz();
tSrcVoxel += srcSampler.peekVoxel1px0py1pz();
tSrcVoxel += srcSampler.peekVoxel1px1py1nz();
tSrcVoxel += srcSampler.peekVoxel1px1py0pz();
tSrcVoxel += srcSampler.peekVoxel1px1py1pz();
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1px1ny1nz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1px1ny0pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1px1ny1pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1px0py1nz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1px0py0pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1px0py1pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1px1py1nz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1px1py0pz());
tSrcVoxel += static_cast<AccumulationType>(srcSampler.peekVoxel1px1py1pz());
tSrcVoxel /= 27;
//tSrcVoxel.setDensity(uDensity);
m_pVolDst->setVoxelAt(iSrcX, iSrcY, iSrcZ, tSrcVoxel);
m_pVolDst->setVoxelAt(iSrcX, iSrcY, iSrcZ, static_cast<DstVolumeType::VoxelType>(tSrcVoxel));
}
}
}
}
template< typename SrcVolumeType, typename DstVolumeType>
void LowPassFilter<SrcVolumeType, DstVolumeType>::executeSAT()
template< typename SrcVolumeType, typename DstVolumeType, typename AccumulationType>
void LowPassFilter<SrcVolumeType, DstVolumeType, AccumulationType>::executeSAT()
{
const uint32_t border = (m_uKernelSize - 1) / 2;
@ -124,7 +122,7 @@ namespace PolyVox
//Use floats for the SAT volume to ensure it works with negative
//densities and with both integral and floating point input volumes.
RawVolume<float> satVolume(Region(satLowerCorner, satUpperCorner));
RawVolume<AccumulationType> satVolume(Region(satLowerCorner, satUpperCorner));
//Clear to zeros (necessary?)
//FIXME - use Volume::fill() method. Implemented in base class as below
@ -140,9 +138,9 @@ namespace PolyVox
}
}
RawVolume<float>::Sampler satVolumeIter(&satVolume);
RawVolume<AccumulationType>::Sampler satVolumeIter(&satVolume);
IteratorController<RawVolume<float>::Sampler> satIterCont;
IteratorController<RawVolume<AccumulationType>::Sampler> satIterCont;
satIterCont.m_regValid = Region(satLowerCorner, satUpperCorner);
satIterCont.m_Iter = &satVolumeIter;
satIterCont.reset();
@ -154,32 +152,34 @@ namespace PolyVox
srcIterCont.m_Iter = &srcVolumeIter;
srcIterCont.reset();
do
/*do
{
float previousSum = satVolumeIter.peekVoxel1nx0py0pz();
AccumulationType previousSum = satVolumeIter.peekVoxel1nx0py0pz();
float currentVal = static_cast<float>(srcVolumeIter.getVoxel().getDensity());
AccumulationType currentVal = srcVolumeIter.getVoxel();
satVolumeIter.setVoxel(previousSum + currentVal);
srcIterCont.moveForward();
}while(satIterCont.moveForward());
}while(satIterCont.moveForward());*/
//Build SAT in three passes
/*for(int32_t z = satLowerCorner.getZ(); z <= satUpperCorner.getZ(); z++)
for(int32_t z = satLowerCorner.getZ(); z <= satUpperCorner.getZ(); z++)
{
for(int32_t y = satLowerCorner.getY(); y <= satUpperCorner.getY(); y++)
{
for(int32_t x = satLowerCorner.getX(); x <= satUpperCorner.getX(); x++)
{
uint32_t previousSum = satVolume.getVoxelAt(x-1,y,z);
uint32_t currentVal = m_pVolSrc->getVoxelAt(x,y,z).getDensity();
AccumulationType previousSum = static_cast<AccumulationType>(satVolume.getVoxelAt(x-1,y,z));
AccumulationType currentVal = static_cast<AccumulationType>(m_pVolSrc->getVoxelAt(x,y,z));
satVolume.setVoxelAt(x,y,z,previousSum + currentVal);
AccumulationType sum = previousSum + currentVal;
satVolume.setVoxelAt(x,y,z,sum);
}
}
}*/
}
for(int32_t z = satLowerCorner.getZ(); z <= satUpperCorner.getZ(); z++)
{
@ -187,8 +187,8 @@ namespace PolyVox
{
for(int32_t x = satLowerCorner.getX(); x <= satUpperCorner.getX(); x++)
{
float previousSum = satVolume.getVoxelAt(x,y-1,z);
float currentSum = satVolume.getVoxelAt(x,y,z);
AccumulationType previousSum = static_cast<AccumulationType>(satVolume.getVoxelAt(x,y-1,z));
AccumulationType currentSum = static_cast<AccumulationType>(satVolume.getVoxelAt(x,y,z));
satVolume.setVoxelAt(x,y,z,previousSum + currentSum);
}
@ -201,8 +201,8 @@ namespace PolyVox
{
for(int32_t x = satLowerCorner.getX(); x <= satUpperCorner.getX(); x++)
{
float previousSum = satVolume.getVoxelAt(x,y,z-1);
float currentSum = satVolume.getVoxelAt(x,y,z);
AccumulationType previousSum = static_cast<AccumulationType>(satVolume.getVoxelAt(x,y,z-1));
AccumulationType currentSum = static_cast<AccumulationType>(satVolume.getVoxelAt(x,y,z));
satVolume.setVoxelAt(x,y,z,previousSum + currentSum);
}
@ -229,27 +229,27 @@ namespace PolyVox
int32_t satUpperY = iSrcY + border;
int32_t satUpperZ = iSrcZ + border;
float a = satVolume.getVoxelAt(satLowerX,satLowerY,satLowerZ);
float b = satVolume.getVoxelAt(satUpperX,satLowerY,satLowerZ);
float c = satVolume.getVoxelAt(satLowerX,satUpperY,satLowerZ);
float d = satVolume.getVoxelAt(satUpperX,satUpperY,satLowerZ);
float e = satVolume.getVoxelAt(satLowerX,satLowerY,satUpperZ);
float f = satVolume.getVoxelAt(satUpperX,satLowerY,satUpperZ);
float g = satVolume.getVoxelAt(satLowerX,satUpperY,satUpperZ);
float h = satVolume.getVoxelAt(satUpperX,satUpperY,satUpperZ);
AccumulationType a = satVolume.getVoxelAt(satLowerX,satLowerY,satLowerZ);
AccumulationType b = satVolume.getVoxelAt(satUpperX,satLowerY,satLowerZ);
AccumulationType c = satVolume.getVoxelAt(satLowerX,satUpperY,satLowerZ);
AccumulationType d = satVolume.getVoxelAt(satUpperX,satUpperY,satLowerZ);
AccumulationType e = satVolume.getVoxelAt(satLowerX,satLowerY,satUpperZ);
AccumulationType f = satVolume.getVoxelAt(satUpperX,satLowerY,satUpperZ);
AccumulationType g = satVolume.getVoxelAt(satLowerX,satUpperY,satUpperZ);
AccumulationType h = satVolume.getVoxelAt(satUpperX,satUpperY,satUpperZ);
float sum = h+c-d-g-f-a+b+e;
AccumulationType sum = h+c-d-g-f-a+b+e;
uint32_t sideLength = border * 2 + 1;
float average = sum / (static_cast<float>(sideLength*sideLength*sideLength));
AccumulationType average = sum / (sideLength*sideLength*sideLength);
//Note: These lines need consideration if src and dest have different voxel types.
typename SrcVolumeType::VoxelType voxel = m_pVolSrc->getVoxelAt(iDstX, iDstY, iDstZ);
//typename SrcVolumeType::VoxelType voxel = m_pVolSrc->getVoxelAt(iDstX, iDstY, iDstZ);
voxel.setDensity(static_cast<typename SrcVolumeType::VoxelType::DensityType>(average));
//voxel.setDensity(static_cast<typename SrcVolumeType::VoxelType::DensityType>(average));
m_pVolDst->setVoxelAt(iDstX, iDstY, iDstZ, voxel);
m_pVolDst->setVoxelAt(iDstX, iDstY, iDstZ, static_cast<DstVolumeType::VoxelType>(average));
//float maxSolid = border * 2/* + 1*/;

View File

@ -58,29 +58,30 @@ void TestLowPassFilter::testExecute()
RawVolume<Density8> resultVolume(reg);
LowPassFilter< RawVolume<Density8>, RawVolume<Density8> > pass1(&volData, reg, &resultVolume, reg, 5);
LowPassFilter< RawVolume<Density8>, RawVolume<Density8>, Density16 > lowPassfilter(&volData, reg, &resultVolume, reg, 5);
pass1.execute();
//Test the normal implementation
lowPassfilter.execute();
QCOMPARE(resultVolume.getVoxelAt(0,0,0), Density8(4));
QCOMPARE(resultVolume.getVoxelAt(1,1,1), Density8(21));
QCOMPARE(resultVolume.getVoxelAt(2,2,2), Density8(10));
QCOMPARE(resultVolume.getVoxelAt(3,3,3), Density8(21));
QCOMPARE(resultVolume.getVoxelAt(4,4,4), Density8(10));
QCOMPARE(resultVolume.getVoxelAt(5,5,5), Density8(21));
QCOMPARE(resultVolume.getVoxelAt(6,6,6), Density8(10));
QCOMPARE(resultVolume.getVoxelAt(7,7,7), Density8(4));
std::cout << "Input volume:" << std::endl;
std::cout << "Voxel = " << static_cast<int>(volData.getVoxelAt(0,0,0).getDensity()) << std::endl; // 32
std::cout << "Voxel = " << static_cast<int>(volData.getVoxelAt(1,1,1).getDensity()) << std::endl; // 0
std::cout << "Voxel = " << static_cast<int>(volData.getVoxelAt(2,2,2).getDensity()) << std::endl; // 3
std::cout << "Voxel = " << static_cast<int>(volData.getVoxelAt(3,3,3).getDensity()) << std::endl; // 0
std::cout << "Voxel = " << static_cast<int>(volData.getVoxelAt(4,4,4).getDensity()) << std::endl; // 32
std::cout << "Voxel = " << static_cast<int>(volData.getVoxelAt(5,5,5).getDensity()) << std::endl; // 0
std::cout << "Voxel = " << static_cast<int>(volData.getVoxelAt(6,6,6).getDensity()) << std::endl; // 32
std::cout << "Voxel = " << static_cast<int>(volData.getVoxelAt(7,7,7).getDensity()) << std::endl; // 0
std::cout << std::endl << "Output volume:" << std::endl;
std::cout << "Voxel = " << static_cast<int>(resultVolume.getVoxelAt(0,0,0).getDensity()) << std::endl; // 4
std::cout << "Voxel = " << static_cast<int>(resultVolume.getVoxelAt(1,1,1).getDensity()) << std::endl; // 21
std::cout << "Voxel = " << static_cast<int>(resultVolume.getVoxelAt(2,2,2).getDensity()) << std::endl; // 10
std::cout << "Voxel = " << static_cast<int>(resultVolume.getVoxelAt(3,3,3).getDensity()) << std::endl; // 21
std::cout << "Voxel = " << static_cast<int>(resultVolume.getVoxelAt(4,4,4).getDensity()) << std::endl; // 10
std::cout << "Voxel = " << static_cast<int>(resultVolume.getVoxelAt(5,5,5).getDensity()) << std::endl; // 21
std::cout << "Voxel = " << static_cast<int>(resultVolume.getVoxelAt(6,6,6).getDensity()) << std::endl; // 10
std::cout << "Voxel = " << static_cast<int>(resultVolume.getVoxelAt(7,7,7).getDensity()) << std::endl; // 4
//Test the SAT implmentation
//FIXME - Shouldn't the results be the same as the normal case?
lowPassfilter.executeSAT();
QCOMPARE(resultVolume.getVoxelAt(0,0,0), Density8(4));
QCOMPARE(resultVolume.getVoxelAt(1,1,1), Density8(8));
QCOMPARE(resultVolume.getVoxelAt(2,2,2), Density8(19));
QCOMPARE(resultVolume.getVoxelAt(3,3,3), Density8(12));
QCOMPARE(resultVolume.getVoxelAt(4,4,4), Density8(19));
QCOMPARE(resultVolume.getVoxelAt(5,5,5), Density8(12));
QCOMPARE(resultVolume.getVoxelAt(6,6,6), Density8(8));
QCOMPARE(resultVolume.getVoxelAt(7,7,7), Density8(2));
}
QTEST_MAIN(TestLowPassFilter)