ITK  6.0.0
Insight Toolkit
SphinxExamples/src/Filtering/Convolution/NormalizedCorrelation/Code.cxx
/*=========================================================================
*
* Copyright NumFOCUS
*
* Licensed under the Apache License, Version 2.0 (the "License");
* you may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
* https://www.apache.org/licenses/LICENSE-2.0.txt
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*
*=========================================================================*/
#include "itkImage.h"
#ifdef ENABLE_QUICKVIEW
# include "QuickView.h"
#endif
#include <iostream>
#include <string>
using FloatImageType = itk::Image<float, 2>;
using UnsignedCharImageType = itk::Image<unsigned char, 2>;
int
main(int argc, char * argv[])
{
if (argc < 2)
{
std::cerr << "Required: filename" << std::endl;
return EXIT_FAILURE;
}
const auto input = itk::ReadImage<FloatImageType>(argv[1]);
// Extract a small region
auto extractFilter = ExtractFilterType::New();
start.Fill(50);
patchSize.Fill(51);
FloatImageType::RegionType desiredRegion(start, patchSize);
extractFilter->SetRegionOfInterest(desiredRegion);
extractFilter->SetInput(input);
extractFilter->Update();
// Perform normalized correlation
// <input type, mask type (not used), output type>
kernelOperator.SetImageKernel(extractFilter->GetOutput());
// The radius of the kernel must be the radius of the patch, NOT the size of the patch
itk::Size<2> radius = extractFilter->GetOutput()->GetLargestPossibleRegion().GetSize();
radius[0] = (radius[0] - 1) / 2;
radius[1] = (radius[1] - 1) / 2;
kernelOperator.CreateToRadius(radius);
auto correlationFilter = CorrelationFilterType::New();
correlationFilter->SetInput(input);
correlationFilter->SetTemplate(kernelOperator);
correlationFilter->Update();
using MinimumMaximumImageCalculatorType = itk::MinimumMaximumImageCalculator<FloatImageType>;
MinimumMaximumImageCalculatorType::Pointer minimumMaximumImageCalculatorFilter =
minimumMaximumImageCalculatorFilter->SetImage(correlationFilter->GetOutput());
minimumMaximumImageCalculatorFilter->Compute();
itk::Index<2> maximumCorrelationPatchCenter = minimumMaximumImageCalculatorFilter->GetIndexOfMaximum();
std::cout << "Maximum: " << maximumCorrelationPatchCenter << std::endl;
// Note that the best correlation is at the center of the patch we extracted (ie. (75, 75) rather than the corner
// (50,50)
{
auto rescaleFilter = RescaleFilterType::New();
rescaleFilter->SetInput(correlationFilter->GetOutput());
rescaleFilter->SetOutputMinimum(0);
rescaleFilter->SetOutputMaximum(255);
rescaleFilter->Update();
itk::WriteImage(rescaleFilter->GetOutput(), "correlation.png");
}
{
auto rescaleFilter = RescaleFilterType::New();
rescaleFilter->SetInput(extractFilter->GetOutput());
rescaleFilter->SetOutputMinimum(0);
rescaleFilter->SetOutputMaximum(255);
rescaleFilter->Update();
itk::WriteImage(rescaleFilter->GetOutput(), "patch.png");
}
// Extract the best matching patch
FloatImageType::IndexType bestPatchStart;
bestPatchStart[0] = maximumCorrelationPatchCenter[0] - radius[0];
bestPatchStart[1] = maximumCorrelationPatchCenter[1] - radius[1];
FloatImageType::RegionType bestPatchRegion(bestPatchStart, patchSize);
auto bestPatchExtractFilter = ExtractFilterType::New();
bestPatchExtractFilter->SetRegionOfInterest(bestPatchRegion);
bestPatchExtractFilter->SetInput(input);
bestPatchExtractFilter->Update();
#ifdef ENABLE_QUICKVIEW
QuickView viewer;
viewer.AddImage(input);
viewer.AddImage(extractFilter->GetOutput());
viewer.AddImage(correlationFilter->GetOutput());
viewer.AddImage(bestPatchExtractFilter->GetOutput());
viewer.Visualize();
#endif
return EXIT_SUCCESS;
}
Pointer
SmartPointer< Self > Pointer
Definition: itkAddImageFilter.h:93
itk::ImageKernelOperator
A NeighborhoodOperator whose coefficients are from an image.
Definition: itkImageKernelOperator.h:49
itk::Index
Represent a n-dimensional index in a n-dimensional image.
Definition: itkIndex.h:68
itk::Size< 2 >
itk::MinimumMaximumImageCalculator
Computes the minimum and the maximum intensity values of an image.
Definition: itkMinimumMaximumImageCalculator.h:45
itkImageFileReader.h
itk::GTest::TypedefsAndConstructors::Dimension2::SizeType
ImageBaseType::SizeType SizeType
Definition: itkGTestTypedefsAndConstructors.h:49
itkImage.h
itkMinimumMaximumImageCalculator.h
itk::Index::Fill
void Fill(IndexValueType value)
Definition: itkIndex.h:272
itk::Size::Fill
void Fill(SizeValueType value)
Definition: itkSize.h:211
itk::GTest::TypedefsAndConstructors::Dimension2::IndexType
ImageBaseType::IndexType IndexType
Definition: itkGTestTypedefsAndConstructors.h:50
QuickView.h
QuickView::AddImage
void AddImage(TImage *, bool FlipVertical=true, std::string Description="")
itk::GTest::TypedefsAndConstructors::Dimension2::RegionType
ImageBaseType::RegionType RegionType
Definition: itkGTestTypedefsAndConstructors.h:54
itkRescaleIntensityImageFilter.h
itkImageFileWriter.h
itk::ImageKernelOperator::SetImageKernel
void SetImageKernel(const ImageType *kernel)
itkRegionOfInterestImageFilter.h
itk::NormalizedCorrelationImageFilter
Computes the normalized correlation of an image and a template.
Definition: itkNormalizedCorrelationImageFilter.h:56
QuickView
A convenient class to render itk images with vtk.
Definition: QuickView.h:111
itk::RescaleIntensityImageFilter
Applies a linear transformation to the intensity levels of the input Image.
Definition: itkRescaleIntensityImageFilter.h:133
itkImageKernelOperator.h
itk::Image
Templated n-dimensional image class.
Definition: itkImage.h:88
New
static Pointer New()
itk::NeighborhoodOperator::CreateToRadius
virtual void CreateToRadius(const SizeType &)
QuickView::Visualize
void Visualize(bool interact=true)
itk::RegionOfInterestImageFilter
Extract a region of interest from the input image.
Definition: itkRegionOfInterestImageFilter.h:54
itkNormalizedCorrelationImageFilter.h
itk::WriteImage
ITK_TEMPLATE_EXPORT void WriteImage(TImagePointer &&image, const std::string &filename, bool compress=false)
Definition: itkImageFileWriter.h:256
itk::Size::GetSize
const SizeValueType * GetSize() const
Definition: itkSize.h:169