Skip to content
Snippets Groups Projects
Commit bc807bce authored by Guillaume Pasero's avatar Guillaume Pasero
Browse files

ENH: implement sub-sampled computation in Advanced and HigherOrder textures

parent 6ccac1f9
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Modules/Feature/Textures/include/otbScalarImageToAdvancedTexturesFilter.h
+ 19
− 0
View file @ bc807bce
......@@ -164,6 +164,18 @@ public:
/** Get the input image maximum */
itkGetMacro(InputImageMaximum, InputPixelType);
/** Set the sub-sampling factor */
itkSetMacro(SubsampleFactor, SizeType);
/** Get the sub-sampling factor */
itkGetMacro(SubsampleFactor, SizeType);
/** Set the sub-sampling offset */
itkSetMacro(SubsampleOffset, OffsetType);
/** Get the sub-sampling offset */
itkGetMacro(SubsampleOffset, OffsetType);
/** Get the mean output image */
OutputImageType * GetMeanOutput();
......@@ -199,6 +211,8 @@ protected:
ScalarImageToAdvancedTexturesFilter();
/** Destructor */
~ScalarImageToAdvancedTexturesFilter() ITK_OVERRIDE;
/** Generate the output informations */
void GenerateOutputInformation() ITK_OVERRIDE;
/** Generate the input requested region */
void GenerateInputRequestedRegion() ITK_OVERRIDE;
/** Before Parallel textures extraction */
......@@ -231,6 +245,11 @@ private:
/** Input image maximum */
InputPixelType m_InputImageMaximum;
/** Sub-sampling factor */
SizeType m_SubsampleFactor;
/** Sub-sampling offset */
OffsetType m_SubsampleOffset;
};
} // End namespace otb
......
Modules/Feature/Textures/include/otbScalarImageToAdvancedTexturesFilter.txx
+ 52
− 1
View file @ bc807bce
......@@ -37,6 +37,8 @@ ScalarImageToAdvancedTexturesFilter<TInputImage, TOutputImage>
, m_NumberOfBinsPerAxis(8)
, m_InputImageMinimum(0)
, m_InputImageMaximum(255)
, m_SubsampleFactor()
, m_SubsampleOffset()
{
// There are 10 outputs corresponding to the 9 textures indices
this->SetNumberOfRequiredOutputs(10);
......@@ -52,6 +54,9 @@ ScalarImageToAdvancedTexturesFilter<TInputImage, TOutputImage>
this->SetNthOutput(7, OutputImageType::New());
this->SetNthOutput(8, OutputImageType::New());
this->SetNthOutput(9, OutputImageType::New());
this->m_SubsampleFactor.Fill(1);
this->m_SubsampleOffset.Fill(0);
}
template <class TInputImage, class TOutputImage>
......@@ -188,6 +193,38 @@ ScalarImageToAdvancedTexturesFilter<TInputImage, TOutputImage>
return static_cast<OutputImageType *>(this->GetOutput(9));
}
template <class TInputImage, class TOutputImage>
void
ScalarImageToAdvancedTexturesFilter<TInputImage, TOutputImage>
::GenerateOutputInformation()
{
// First, call superclass implementation
Superclass::GenerateOutputInformation();
// Compute output size, origin & spacing
InputRegionType inputRegion = this->GetInput()->GetLargestPossibleRegion();
OutputRegionType outputRegion;
outputRegion.SetIndex(0,0);
outputRegion.SetIndex(1,0);
outputRegion.SetSize(0, 1 + (inputRegion.GetSize(0) - 1 - m_SubsampleOffset[0]) / m_SubsampleFactor[0]);
outputRegion.SetSize(1, 1 + (inputRegion.GetSize(1) - 1 - m_SubsampleOffset[1]) / m_SubsampleFactor[1]);
typename OutputImageType::SpacingType outSpacing = this->GetInput()->GetSpacing();
outSpacing[0] *= m_SubsampleFactor[0];
outSpacing[1] *= m_SubsampleFactor[1];
typename OutputImageType::PointType outOrigin;
this->GetInput()->TransformIndexToPhysicalPoint(inputRegion.GetIndex()+m_SubsampleOffset,outOrigin);
for (unsigned int i=0 ; i<this->GetNumberOfOutputs() ; i++)
{
OutputImagePointerType outputPtr = this->GetOutput(i);
outputPtr->SetLargestPossibleRegion(outputRegion);
outputPtr->SetOrigin(outOrigin);
outputPtr->SetSpacing(outSpacing);
}
}
template <class TInputImage, class TOutputImage>
void
ScalarImageToAdvancedTexturesFilter<TInputImage, TOutputImage>
......@@ -209,12 +246,19 @@ ScalarImageToAdvancedTexturesFilter<TInputImage, TOutputImage>
// We use only the first output since requested regions for all outputs are enforced to be equal
// by the default GenerateOutputRequestedRegiont() implementation
OutputRegionType outputRequestedRegion = outputPtr->GetRequestedRegion();
InputRegionType inputLargest = inputPtr->GetLargestPossibleRegion();
typename OutputRegionType::IndexType outputIndex = outputRequestedRegion.GetIndex();
typename OutputRegionType::SizeType outputSize = outputRequestedRegion.GetSize();
typename InputRegionType::IndexType inputIndex;
typename InputRegionType::SizeType inputSize;
// Convert index and size to full grid
outputIndex[0] = outputIndex[0] * m_SubsampleFactor[0] + m_SubsampleOffset[0] + inputLargest.GetIndex(0);
outputIndex[1] = outputIndex[1] * m_SubsampleFactor[1] + m_SubsampleOffset[1] + inputLargest.GetIndex(1);
outputSize[0] = 1 + (outputSize[0] - 1) * m_SubsampleFactor[0];
outputSize[1] = 1 + (outputSize[1] - 1) * m_SubsampleFactor[1];
// First, apply offset
for (unsigned int dim = 0; dim < InputImageType::ImageDimension; ++dim)
{
......@@ -311,6 +355,8 @@ ScalarImageToAdvancedTexturesFilter<TInputImage, TOutputImage>
const unsigned int histSize = m_NumberOfBinsPerAxis;
const long unsigned int twiceHistSize = 2 * m_NumberOfBinsPerAxis;
InputRegionType inputLargest = inputPtr->GetLargestPossibleRegion();
// Set-up progress reporting
itk::ProgressReporter progress(this, threadId, outputRegionForThread.GetNumberOfPixels());
......@@ -330,11 +376,16 @@ ScalarImageToAdvancedTexturesFilter<TInputImage, TOutputImage>
typename InputRegionType::IndexType inputIndex;
typename InputRegionType::SizeType inputSize;
// Convert index to full grid
typename OutputImageType::IndexType outIndex;
// First, create an window for neighborhood iterator based on m_Radius
// For example, if xradius and yradius is 2. window size is 5x5 (2 * radius + 1).
for (unsigned int dim = 0; dim < InputImageType::ImageDimension; ++dim)
{
inputIndex[dim] = varianceIt.GetIndex()[dim] - m_Radius[dim];
outIndex[dim] = varianceIt.GetIndex()[dim] * m_SubsampleFactor[dim]
+ m_SubsampleOffset[dim] + inputLargest.GetIndex(dim);
inputIndex[dim] = outIndex[dim] - m_Radius[dim];
inputSize[dim] = 2 * m_Radius[dim] + 1;
}
......
Modules/Feature/Textures/include/otbScalarImageToHigherOrderTexturesFilter.h
+ 20
− 0
View file @ bc807bce
......@@ -148,6 +148,18 @@ public:
itkSetMacro(FastCalculations, bool);
itkBooleanMacro(FastCalculations);
/** Set the sub-sampling factor */
itkSetMacro(SubsampleFactor, SizeType);
/** Get the sub-sampling factor */
itkGetMacro(SubsampleFactor, SizeType);
/** Set the sub-sampling offset */
itkSetMacro(SubsampleOffset, OffsetType);
/** Get the sub-sampling offset */
itkGetMacro(SubsampleOffset, OffsetType);
/** Get the Short Run Emphasis output image */
OutputImageType * GetShortRunEmphasisOutput();
......@@ -186,6 +198,8 @@ protected:
ScalarImageToHigherOrderTexturesFilter();
/** Destructor */
~ScalarImageToHigherOrderTexturesFilter() ITK_OVERRIDE;
/** Generate the output informations */
void GenerateOutputInformation() ITK_OVERRIDE;
/** Generate the input requested region */
void GenerateInputRequestedRegion() ITK_OVERRIDE;
/** Parallel textures extraction */
......@@ -215,6 +229,12 @@ private:
/** Fast calculation */
bool m_FastCalculations;
/** Sub-sampling factor */
SizeType m_SubsampleFactor;
/** Sub-sampling offset */
OffsetType m_SubsampleOffset;
};
} // End namespace otb
......
Modules/Feature/Textures/include/otbScalarImageToHigherOrderTexturesFilter.txx
+ 57
− 4
View file @ bc807bce
......@@ -32,7 +32,9 @@ ScalarImageToHigherOrderTexturesFilter<TInputImage, TOutputImage>
m_NumberOfBinsPerAxis(8),
m_InputImageMinimum(0),
m_InputImageMaximum(255),
m_FastCalculations(false)
m_FastCalculations(false),
m_SubsampleFactor(),
m_SubsampleOffset()
{
// There are 11 outputs corresponding to the 8 textures indices
this->SetNumberOfRequiredOutputs(10);
......@@ -71,6 +73,8 @@ ScalarImageToHigherOrderTexturesFilter<TInputImage, TOutputImage>
}
this->SetOffsets( offsets );
this->m_SubsampleFactor.Fill(1);
this->m_SubsampleOffset.Fill(0);
}
template <class TInputImage, class TOutputImage>
......@@ -231,6 +235,38 @@ ScalarImageToHigherOrderTexturesFilter<TInputImage, TOutputImage>
this->SetOffsets( offsetVector );
}
template <class TInputImage, class TOutputImage>
void
ScalarImageToHigherOrderTexturesFilter<TInputImage, TOutputImage>
::GenerateOutputInformation()
{
// First, call superclass implementation
Superclass::GenerateOutputInformation();
// Compute output size, origin & spacing
InputRegionType inputRegion = this->GetInput()->GetLargestPossibleRegion();
OutputRegionType outputRegion;
outputRegion.SetIndex(0,0);
outputRegion.SetIndex(1,0);
outputRegion.SetSize(0, 1 + (inputRegion.GetSize(0) - 1 - m_SubsampleOffset[0]) / m_SubsampleFactor[0]);
outputRegion.SetSize(1, 1 + (inputRegion.GetSize(1) - 1 - m_SubsampleOffset[1]) / m_SubsampleFactor[1]);
typename OutputImageType::SpacingType outSpacing = this->GetInput()->GetSpacing();
outSpacing[0] *= m_SubsampleFactor[0];
outSpacing[1] *= m_SubsampleFactor[1];
typename OutputImageType::PointType outOrigin;
this->GetInput()->TransformIndexToPhysicalPoint(inputRegion.GetIndex()+m_SubsampleOffset,outOrigin);
for (unsigned int i=0 ; i<this->GetNumberOfOutputs() ; i++)
{
OutputImagePointerType outputPtr = this->GetOutput(i);
outputPtr->SetLargestPossibleRegion(outputRegion);
outputPtr->SetOrigin(outOrigin);
outputPtr->SetSpacing(outSpacing);
}
}
template <class TInputImage, class TOutputImage>
void
ScalarImageToHigherOrderTexturesFilter<TInputImage, TOutputImage>
......@@ -252,7 +288,17 @@ ScalarImageToHigherOrderTexturesFilter<TInputImage, TOutputImage>
// We use only the first output since requested regions for all outputs are enforced to be equal
// by the default GenerateOutputRequestedRegiont() implementation
OutputRegionType outputRequestedRegion = outputPtr->GetRequestedRegion();
InputRegionType inputRequestedRegion = outputRequestedRegion;
typename OutputRegionType::IndexType outputIndex = outputRequestedRegion.GetIndex();
typename OutputRegionType::SizeType outputSize = outputRequestedRegion.GetSize();
InputRegionType inputLargest = inputPtr->GetLargestPossibleRegion();
// Convert index and size to full grid
outputIndex[0] = outputIndex[0] * m_SubsampleFactor[0] + m_SubsampleOffset[0] + inputLargest.GetIndex(0);
outputIndex[1] = outputIndex[1] * m_SubsampleFactor[1] + m_SubsampleOffset[1] + inputLargest.GetIndex(1);
outputSize[0] = 1 + (outputSize[0] - 1) * m_SubsampleFactor[0];
outputSize[1] = 1 + (outputSize[1] - 1) * m_SubsampleFactor[1];
InputRegionType inputRequestedRegion(outputIndex,outputSize);
// Apply the radius
inputRequestedRegion.PadByRadius(m_Radius);
......@@ -297,6 +343,8 @@ ScalarImageToHigherOrderTexturesFilter<TInputImage, TOutputImage>
inputPtr->TransformIndexToPhysicalPoint( outputImagesIterators[0].GetIndex() + m_Radius, bottomRightPoint );
double maxDistance = topLeftPoint.EuclideanDistanceTo(bottomRightPoint);
InputRegionType inputLargest = inputPtr->GetLargestPossibleRegion();
// Set-up progress reporting
itk::ProgressReporter progress(this, threadId, outputRegionForThread.GetNumberOfPixels());
......@@ -304,12 +352,17 @@ ScalarImageToHigherOrderTexturesFilter<TInputImage, TOutputImage>
while ( !outputImagesIterators[0].IsAtEnd() )
{
// Compute the region on which run-length matrix will be estimated
typename InputRegionType::IndexType inputIndex = outputImagesIterators[0].GetIndex() - m_Radius;
typename InputRegionType::IndexType inputIndex;
typename InputRegionType::SizeType inputSize;
// First, apply offset
// Convert index to full grid
typename OutputImageType::IndexType outIndex;
for (unsigned int dim = 0; dim < InputImageType::ImageDimension; ++dim)
{
outIndex[dim] = outputImagesIterators[0].GetIndex()[dim] * m_SubsampleFactor[dim]
+ m_SubsampleOffset[dim] + inputLargest.GetIndex(dim);
inputIndex[dim] = outIndex[dim] - m_Radius[dim];
inputSize[dim] = 2 * m_Radius[dim] + 1;
}
......
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