Commit 525a676b authored by Julien Malik's avatar Julien Malik

ENH: add ImageRegionSquareTileSplitter

parent 6d5d9554
......@@ -104,29 +104,27 @@ public:
/** Region typedef support. */
typedef itk::ImageRegion<VImageDimension> RegionType;
itkSetMacro(PixelSizeInBytes, unsigned int);
itkGetMacro(PixelSizeInBytes, unsigned int);
itkSetMacro(TileSizeInBytes, unsigned int);
itkGetMacro(TileSizeInBytes, unsigned int);
/** How many pieces can the specifed region be split? A given region
* cannot always be divided into the requested number of pieces. For
* instance, if the numberOfPieces exceeds the number of pixels along
* a certain dimensions, then some splits will not be possible. This
* method returns a number less than or equal to the requested number
* of pieces. */
/** How many pieces can the specified region be split? A given region
* cannot always be divided into the requested number of pieces. For
* instance, if the numberOfPieces exceeds the number of pixels along
* a certain dimensions, then some splits will not be possible.
*/
virtual unsigned int GetNumberOfSplits(const RegionType& region,
unsigned int requestedNumber);
/** Get a region definition that represents the ith piece a specified region.
* The "numberOfPieces" specified should be less than or equal to what
* The "numberOfPieces" specified should be less than or equal to what
* GetNumberOfSplits() returns. */
virtual RegionType GetSplit(unsigned int i, unsigned int numberOfPieces,
const RegionType& region);
itkGetMacro(TileSizeAlignment, unsigned int);
itkSetMacro(TileSizeAlignment, unsigned int);
itkGetMacro(TileDimension, unsigned int);
protected:
ImageRegionSquareTileSplitter() : m_AlignStep(0), m_PixelSizeInBytes(0), m_TileSizeInBytes(0){}
ImageRegionSquareTileSplitter() : m_SplitsPerDimension(0U), m_TileDimension(0), m_TileSizeAlignment(16) {}
virtual ~ImageRegionSquareTileSplitter() {}
void PrintSelf(std::ostream& os, itk::Indent indent) const;
......@@ -134,10 +132,9 @@ private:
ImageRegionSquareTileSplitter(const ImageRegionSquareTileSplitter &); //purposely not implemented
void operator =(const ImageRegionSquareTileSplitter&); //purposely not implemented
unsigned int m_SplitsPerDimension[VImageDimension];
unsigned int m_AlignStep;
unsigned int m_PixelSizeInBytes;
unsigned int m_TileSizeInBytes;
itk::FixedArray<unsigned int, VImageDimension> m_SplitsPerDimension;
unsigned int m_TileDimension;
unsigned int m_TileSizeAlignment;
};
} // end namespace otb
......
......@@ -25,140 +25,77 @@
namespace otb
{
/**
*
*/
template <unsigned int VImageDimension>
unsigned int
ImageRegionSquareTileSplitter<VImageDimension>
::GetNumberOfSplits(const RegionType& region, unsigned int requestedNumber)
{
if (m_TileSizeInBytes == 0)
{
itkExceptionMacro(<< "SetTileSizeInBytes has not been called");
}
if (m_PixelSizeInBytes == 0)
{
itkExceptionMacro(<< "SetPixelSizeInBytes has not been called");
}
unsigned int theoricalNbPixelPerTile = region.GetNumberOfPixels() / requestedNumber;
unsigned int theoricalTileDimension = static_cast<unsigned int> (vcl_sqrt(static_cast<double>(theoricalNbPixelPerTile)) );
// Compute the theorical square tile dimensions
const unsigned int nbPixelPerTile = m_TileSizeInBytes / m_PixelSizeInBytes;
unsigned int nbPixelPerDim = static_cast<unsigned int>( vcl_sqrt(nbPixelPerTile));
// Take the previous multiple of m_TileSizeAlignment (eventually generate more splits than requested)
m_TileDimension = theoricalTileDimension / m_TileSizeAlignment * m_TileSizeAlignment;
// Use a reference tile size corresponding to a 256 * 256 region of a 4 bands unsigned short image (=512 k)
const unsigned int ReferenceTileSizeInBytes = 256 * 256 * 4 * 2;
unsigned int nbPixelPerReferenceTile = ReferenceTileSizeInBytes / m_PixelSizeInBytes;
unsigned int referenceNbPixelPerDim = static_cast<unsigned int>( vcl_sqrt(static_cast<float>(nbPixelPerReferenceTile)) );
// Align the tile dimension to the next multiple of 16 using integer division (TIFF tiles are aligned with 16)
referenceNbPixelPerDim = ( referenceNbPixelPerDim + 15 ) / 16 * 16;
// Align nbPixelPerDim with the next multiple of referenceNbPixelPerDim using integer division
m_AlignStep = (nbPixelPerDim + referenceNbPixelPerDim) / referenceNbPixelPerDim * referenceNbPixelPerDim;
// Use at least 16*16 tiles
if (m_AlignStep < 16)
// Minimal tile size is m_TileSizeAlignment * m_TileSizeAlignment
if (m_TileDimension < m_TileSizeAlignment)
{
m_AlignStep = 16;
otbMsgDevMacro(<< "Warning: clamping tile size to " << m_TileSizeAlignment << " * " << m_TileSizeAlignment);
m_TileDimension = m_TileSizeAlignment;
}
const SizeType& regionSize = region.GetSize();
const IndexType& regionIndex = region.GetIndex();
// requested number of splits per dimension
unsigned int numPieces = 1;
// determine the actual number of pieces that will be generated
for (unsigned int j = VImageDimension; j > 0; --j)
const SizeType& regionSize = region.GetSize();
for (unsigned int j = 0; j < VImageDimension; ++j)
{
// otbMsgDevMacro(<< "*** Dimension: " << j-1);
unsigned long int remainingToDo =
static_cast<unsigned long int>(vcl_ceil(static_cast<double>(requestedNumber) / numPieces));
unsigned int maxPieces = (regionIndex[j - 1] + regionSize[j - 1] - 1) / m_AlignStep - regionIndex[j - 1]
/ m_AlignStep + 1;
unsigned int stepPerPiece = 1;
if (remainingToDo < maxPieces)
{
stepPerPiece = static_cast<unsigned int> (vcl_floor(static_cast<double> (maxPieces) / remainingToDo));
if ((remainingToDo - 1) * (stepPerPiece + 1) < maxPieces)
{
stepPerPiece += 1;
}
}
unsigned int maxPieceUsed = static_cast<unsigned int> (vcl_ceil(static_cast<double> (maxPieces) / stepPerPiece));
m_SplitsPerDimension[j - 1] = maxPieceUsed;
// otbMsgDevMacro("*** maxPieces stepPerPiece maxPieceUsed " << maxPieces
// << " " << stepPerPiece << " " << maxPieceUsed);
numPieces *= maxPieceUsed;
m_SplitsPerDimension[j] = (regionSize[j] + m_TileDimension - 1) / m_TileDimension;
numPieces *= m_SplitsPerDimension[j];
}
// otbMsgDevMacro("*** numPieces " << numPieces);
otbMsgDevMacro(<< "Tile dimension : " << m_TileDimension)
otbMsgDevMacro(<< "Number of splits per dimension : " << m_SplitsPerDimension[0] << " " << m_SplitsPerDimension[1])
return numPieces;
}
/**
*
*/
template <unsigned int VImageDimension>
itk::ImageRegion<VImageDimension>
ImageRegionSquareTileSplitter<VImageDimension>
::GetSplit(unsigned int i, unsigned int numberOfPieces, const RegionType& region)
::GetSplit(unsigned int i, unsigned int itkNotUsed(numberOfPieces), const RegionType& region)
{
RegionType splitRegion;
IndexType splitIndex, regionIndex;
SizeType splitSize, regionSize;
IndexType splitIndex;
// Initialize the splitRegion to the requested region
splitRegion = region;
splitIndex = splitRegion.GetIndex();
splitSize = splitRegion.GetSize();
// Compute the actual number of splits
unsigned int numPieces = 1;
for (unsigned int j = 0; j < VImageDimension; ++j)
{
numPieces *= m_SplitsPerDimension[j];
}
regionSize = region.GetSize();
regionIndex = region.GetIndex();
// Sanity check
if (i >= numPieces)
{
itkExceptionMacro("Asked for split number " << i << " but region contains only " << numPieces << " splits");
}
unsigned int numPieces = GetNumberOfSplits(region, numberOfPieces);
if (i > numPieces)
// Compute the split index in the streaming grid
unsigned int remaining = i;
for (unsigned int j = VImageDimension - 1; j > 0; --j)
{
itkDebugMacro(" Cannot Split");
return splitRegion;
splitIndex[j] = remaining / m_SplitsPerDimension[VImageDimension - 1 - j];
remaining = remaining % m_SplitsPerDimension[VImageDimension - 1 - j];
}
splitIndex[0] = remaining;
unsigned int stackSize = 1;
// Transform the split index to the actual coordinates
for (unsigned int j = 0; j < VImageDimension; ++j)
{
unsigned int slicePos = (i % (stackSize * m_SplitsPerDimension[j])) / stackSize;
stackSize *= m_SplitsPerDimension[j];
unsigned int generalSplitSize =
static_cast<unsigned int> (vcl_ceil(static_cast<double> (regionSize[j]) / (m_SplitsPerDimension[j]
*
m_AlignStep))) * m_AlignStep;
if (slicePos == 0)
{
splitIndex[j] = regionIndex[j];
}
else
{
splitIndex[j] = (regionIndex[j] / generalSplitSize + slicePos) * generalSplitSize;
}
if (slicePos == 0)
{
splitSize[j] = generalSplitSize - (regionIndex[j] % generalSplitSize);
}
else if (slicePos == m_SplitsPerDimension[j] - 1)
{
splitSize[j] = regionSize[j] - (generalSplitSize - (regionIndex[j] % generalSplitSize))
- (m_SplitsPerDimension[j] - 2) * generalSplitSize;
}
else
{
splitSize[j] = generalSplitSize;
}
splitRegion.SetIndex(j, region.GetIndex(j) + m_TileDimension * splitIndex[j]);
splitRegion.SetSize(j, m_TileDimension);
}
// set the split region ivars
splitRegion.SetIndex(splitIndex);
splitRegion.SetSize(splitSize);
// Handle the borders
splitRegion.Crop(region);
return splitRegion;
}
......@@ -172,9 +109,10 @@ ImageRegionSquareTileSplitter<VImageDimension>
::PrintSelf(std::ostream& os, itk::Indent indent) const
{
Superclass::PrintSelf(os, indent);
os << indent << "AlignStep : " << m_AlignStep << std::endl;
os << indent << "PixelSizeInBytes : " << m_PixelSizeInBytes << std::endl;
os << indent << "TileSizeInBytes : " << m_TileSizeInBytes << std::endl;
os << indent << "SplitsPerDimension : " << m_SplitsPerDimension << std::endl;
os << indent << "TileDimension : " << m_TileDimension << std::endl;
os << indent << "TileSizeAlignment : " << m_TileSizeAlignment << std::endl;
}
} // end namespace itk
......
......@@ -743,10 +743,10 @@ ADD_TEST(coTuImageRegionSquareTileSplitterNew ${COMMON_TESTS8}
ADD_TEST(coTvImageRegionSquareTileSplitter ${COMMON_TESTS8}
--compare-ascii ${NOTOL}
${BASELINE_FILES}/coImageRegionTileMapSplitter.txt
${TEMP}/coImageRegionTileMapSplitter.txt
${BASELINE_FILES}/coImageRegionSquareTileSplitter.txt
${TEMP}/coImageRegionSquareTileSplitter.txt
otbImageRegionSquareTileSplitter
${TEMP}/coImageRegionTileMapSplitter.txt
${TEMP}/coImageRegionSquareTileSplitter.txt
)
# ------------- otb::ImageOfVectorsToMonoChannelExtractROI ----------------------------
......
......@@ -19,160 +19,126 @@
#include <fstream>
const int Dimension = 2;
typedef otb::ImageRegionSquareTileSplitter<Dimension> SplitterType;
typedef otb::ImageRegionSquareTileSplitter<Dimension> SquareTileSplitterType;
typedef SquareTileSplitterType::IndexType IndexType;
typedef SquareTileSplitterType::SizeType SizeType;
typedef SquareTileSplitterType::RegionType RegionType;
int otbImageRegionSquareTileSplitterNew(int argc, char * argv[])
{
SplitterType::Pointer splitter = SplitterType::New();
SquareTileSplitterType::Pointer splitter = SquareTileSplitterType::New();
std::cout << splitter << std::endl;
return EXIT_SUCCESS;
}
int otbImageRegionSquareTileSplitter(int argc, char * argv[])
int TestSplitter(const RegionType& region, unsigned int PixelSize, unsigned int MaxTileSize, std::ostream& os)
{
typedef SplitterType::IndexType IndexType;
typedef SplitterType::SizeType SizeType;
typedef SplitterType::RegionType RegionType;
os << "----------------------------------" << std::endl;
os << "Region : " << region << std::endl;
os << "PixelSize : " << PixelSize << std::endl;
os << "MaxTileSize : " << MaxTileSize << std::endl;
SquareTileSplitterType::Pointer splitter;
splitter = SquareTileSplitterType::New();
unsigned int requestedNbSplits = region.GetNumberOfPixels() * PixelSize / MaxTileSize;
if (requestedNbSplits == 0)
requestedNbSplits = 1;
os << "Requested Number of splits : " << requestedNbSplits << std::endl;
const unsigned int nbSplits = splitter->GetNumberOfSplits(region, requestedNbSplits);
os << "Actual Number of splits : " << nbSplits << std::endl;
RegionType split;
// First split :
split = splitter->GetSplit(0, nbSplits, region);
os << "First Split : " << split
<< "(" << split.GetNumberOfPixels() * PixelSize << " bytes)" << std::endl;
if (nbSplits > 1)
{
// Second split :
split = splitter->GetSplit(1, nbSplits, region);
os << "Second Split : " << split
<< "(" << split.GetNumberOfPixels() * PixelSize << " bytes)" << std::endl;
}
if (nbSplits > 2)
{
// Last split :
split = splitter->GetSplit(nbSplits - 1, nbSplits, region);
os << "Last Split : " << split
<< "(" << split.GetNumberOfPixels() * PixelSize << " bytes)" << std::endl;
}
typedef unsigned short PixelPrecisionType;
const unsigned int NbComponents = 8;
const unsigned int PixelSizeInBytes = NbComponents * sizeof(PixelPrecisionType);
return EXIT_SUCCESS;
}
int otbImageRegionSquareTileSplitter(int argc, char * argv[])
{
std::ofstream outfile(argv[1]);
SplitterType::Pointer filter = FilterType::New();
RegionType region, region2;
unsigned int nb, nbSplitTheoric, nbAsked;
IndexType index;
SizeType size;
//Case 1
index[0] = 45;
index[1] = 45;
size[0] = 1000;
size[1] = 1500;
nbSplitTheoric = 10;
nbAsked = 2;
region.SetSize(size);
region.SetIndex(index);
nb = filter->GetNumberOfSplits(region, nbSplitTheoric);
region2 = filter->GetSplit(nbAsked, nb, region);
outfile << "\nCase 1 \n";
outfile << "Input region: " << region << std::endl;
outfile << "Input NumberOfSplits: " << nbSplitTheoric << std::endl;
outfile << "Output GetNumberOfSplits: " << nb << std::endl;
outfile << "Output GetSplit(" << nbAsked << "," << nb << ", input region): " << std::endl;
outfile << "Output region: " << region2 << std::endl;
//Case 2
index[0] = 45;
index[1] = 45;
size[0] = 1048576;
size[1] = 1048576;
nbSplitTheoric = 16777216;
nbAsked = 2;
region.SetSize(size);
region.SetIndex(index);
nb = filter->GetNumberOfSplits(region, nbSplitTheoric);
region2 = filter->GetSplit(nbAsked, nb, region);
outfile << "\nCase 2 \n";
outfile << "Input region: " << region << std::endl;
outfile << "Input NumberOfSplits: " << nbSplitTheoric << std::endl;
outfile << "Output GetNumberOfSplits: " << nb << std::endl;
outfile << "Output GetSplit(" << nbAsked << "," << nb << ", input region): " << std::endl;
outfile << "Output region: " << region2 << std::endl;
//Case 3
index[0] = 45;
index[1] = 45;
size[0] = 1048576;
size[1] = 1048576;
nbSplitTheoric = 23;
nbAsked = 4;
region.SetSize(size);
region.SetIndex(index);
nb = filter->GetNumberOfSplits(region, nbSplitTheoric);
region2 = filter->GetSplit(nbAsked, nb, region);
outfile << "\nCase 3 \n";
outfile << "Input region: " << region << std::endl;
outfile << "Input NumberOfSplits: " << nbSplitTheoric << std::endl;
outfile << "Output GetNumberOfSplits: " << nb << std::endl;
outfile << "Output GetSplit(" << nbAsked << "," << nb << ", input region): " << std::endl;
outfile << "Output region: " << region2 << std::endl;
//Case 4
index[0] = 45;
index[1] = 45;
size[0] = 1048576;
size[1] = 1024;
nbSplitTheoric = 16777216;
nbAsked = 16387;
region.SetSize(size);
region.SetIndex(index);
nb = filter->GetNumberOfSplits(region, nbSplitTheoric);
region2 = filter->GetSplit(nbAsked, nb, region);
outfile << "\nCase 4 \n";
outfile << "Input region: " << region << std::endl;
outfile << "Input NumberOfSplits: " << nbSplitTheoric << std::endl;
outfile << "Output GetNumberOfSplits: " << nb << std::endl;
outfile << "Output GetSplit(" << nbAsked << "," << nb << ", input region): " << std::endl;
outfile << "Output region: " << region2 << std::endl;
//Case 5
index[0] = 0;
index[1] = 0;
size[0] = 513;
size[1] = 5376;
nbSplitTheoric = 8;
nbAsked = 9;
region.SetSize(size);
region.SetIndex(index);
nb = filter->GetNumberOfSplits(region, nbSplitTheoric);
region2 = filter->GetSplit(nbAsked, nb, region);
outfile << "\nCase 5 \n";
outfile << "Input region: " << region << std::endl;
outfile << "Input NumberOfSplits: " << nbSplitTheoric << std::endl;
outfile << "Output GetNumberOfSplits: " << nb << std::endl;
outfile << "Output GetSplit(" << nbAsked << "," << nb << ", input region): " << std::endl;
outfile << "Output region: " << region2 << std::endl;
//Case 6
index[0] = 0;
index[1] = 0;
size[0] = 3;
size[1] = 2;
nbSplitTheoric = 5;
nbAsked = 0;
region.SetSize(size);
region.SetIndex(index);
nb = filter->GetNumberOfSplits(region, nbSplitTheoric);
region2 = filter->GetSplit(nbAsked, nb, region);
outfile << "\nCase 5 \n";
outfile << "Input region: " << region << std::endl;
outfile << "Input NumberOfSplits: " << nbSplitTheoric << std::endl;
outfile << "Output GetNumberOfSplits: " << nb << std::endl;
outfile << "Output GetSplit(" << nbAsked << "," << nb << ", input region): " << std::endl;
outfile << "Output region: " << region2 << std::endl;
RegionType region;
// Test with a 0-based indexed region
region.SetIndex(0, 0);
region.SetIndex(1, 0);
region.SetSize(0, 1024);
region.SetSize(1, 1024);
TestSplitter(region, 1, 128, outfile);
TestSplitter(region, 1, 512*512, outfile);
TestSplitter(region, 2, 512*512, outfile);
TestSplitter(region, 4, 512*512, outfile);
TestSplitter(region, 8, 512*512, outfile);
// Test with a shifted region
region.SetIndex(0, 42);
region.SetIndex(1, 42);
region.SetSize(0, 1000);
region.SetSize(1, 1000);
TestSplitter(region, 1, 128, outfile);
TestSplitter(region, 1, 512*512, outfile);
TestSplitter(region, 2, 512*512, outfile);
TestSplitter(region, 4, 512*512, outfile);
TestSplitter(region, 8, 512*512, outfile);
// Test with a negative shift
region.SetIndex(0, -42);
region.SetIndex(1, -42);
region.SetSize(0, 1000);
region.SetSize(1, 1000);
TestSplitter(region, 1, 128, outfile);
TestSplitter(region, 1, 512*512, outfile);
TestSplitter(region, 2, 512*512, outfile);
TestSplitter(region, 4, 512*512, outfile);
TestSplitter(region, 8, 512*512, outfile);
// Test with a reduced size
region.SetIndex(0, 0);
region.SetIndex(1, 0);
region.SetSize(0, 1);
region.SetSize(1, 1);
TestSplitter(region, 1, 128, outfile);
TestSplitter(region, 1, 512*512, outfile);
TestSplitter(region, 2, 512*512, outfile);
TestSplitter(region, 4, 512*512, outfile);
TestSplitter(region, 8, 512*512, outfile);
// Test with a reduced size, shifted
region.SetIndex(0, 42);
region.SetIndex(1, 42);
region.SetSize(0, 1);
region.SetSize(1, 1);
TestSplitter(region, 1, 128, outfile);
TestSplitter(region, 1, 512*512, outfile);
TestSplitter(region, 2, 512*512, outfile);
TestSplitter(region, 4, 512*512, outfile);
TestSplitter(region, 8, 512*512, outfile);
outfile.close();
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment