HarrisExample.cxx 6.31 KB
 Sébastien Dinot committed Mar 08, 2017 1 /* Julien Osman committed Apr 23, 2020 2 * Copyright (C) 2005-2020 Centre National d'Etudes Spatiales (CNES) Sébastien Dinot committed Mar 08, 2017 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 * * This file is part of Orfeo Toolbox * * https://www.orfeo-toolbox.org/ * * 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 * * http://www.apache.org/licenses/LICENSE-2.0 * * 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. */ Jordi Inglada committed Apr 04, 2006 20 Patrick Imbo committed Jun 29, 2006 21 Guillaume Pasero committed Mar 20, 2019 22 #include "otbMacro.h" Jordi Inglada committed May 18, 2006 23 #include "otbImage.h" Jordi Inglada committed Apr 04, 2006 24 25 26 27 #include "otbImageFileReader.h" #include "otbImageFileWriter.h" Victor Poughon committed Feb 20, 2019 28 29 30 31 /* Example usage: ./HarrisExample Input/ROISpot5.png Output/ROISpot5Harris.png 1.5 2 0.1 */ Jordi Inglada committed Apr 04, 2006 32 Jordi Inglada committed Jul 03, 2006 33 // This example illustrates the use of the \doxygen{otb}{HarrisImageFilter}. Jordi Inglada committed Apr 04, 2006 34 // Emmanuel Christophe committed Dec 06, 2008 35 // The first step required to use this filter is to include its header file. Jordi Inglada committed Apr 04, 2006 36 Jordi Inglada committed Apr 13, 2006 37 #include "otbHarrisImageToPointSetFilter.h" Jordi Inglada committed Apr 04, 2006 38 39 #include "itkRescaleIntensityImageFilter.h" Victor Poughon committed Feb 20, 2019 40 int main(int argc, char* argv[]) Jordi Inglada committed Apr 04, 2006 41 { OTB Bot committed Apr 01, 2010 42 if (argc != 6) Victor Poughon committed Feb 20, 2019 43 { Jordi Inglada committed Apr 04, 2006 44 std::cerr << "Usage: " << argv[0] << " inputImageFile "; Emmanuel Christophe committed Dec 06, 2008 45 std::cerr << " outputHarrisImageFile sigmaD sigmaI alpha" << std::endl; Jordi Inglada committed Apr 04, 2006 46 return EXIT_FAILURE; Victor Poughon committed Feb 20, 2019 47 } Emmanuel Christophe committed Jan 31, 2009 48 Victor Poughon committed Feb 20, 2019 49 50 const char* inputFilename = argv[1]; const char* outputFilename = argv[2]; Emmanuel Christophe committed Jan 31, 2009 51 Victor Poughon committed Feb 20, 2019 52 53 54 double SigmaD((double)::atof(argv[3])); double SigmaI((double)::atof(argv[4])); double Alpha((double)::atof(argv[5])); Emmanuel Christophe committed Jan 31, 2009 55 Victor Poughon committed May 16, 2019 56 57 58 using InputPixelType = float; const unsigned int Dimension = 2; using OutputPixelType = unsigned char; Emmanuel Christophe committed Jan 31, 2009 59 Victor Poughon committed May 16, 2019 60 61 using InputImageType = otb::Image; using OutputImageType = otb::Image; Emmanuel Christophe committed Jan 31, 2009 62 Victor Poughon committed May 16, 2019 63 using ReaderType = otb::ImageFileReader; Emmanuel Christophe committed Jan 31, 2009 64 65 66 67 68 // The \doxygen{otb}{HarrisImageFilter} is templated over the // input and output image types, so we start by // defining: Victor Poughon committed May 16, 2019 69 70 using HarrisFilterType = otb::HarrisImageFilter; using RescalerType = itk::RescaleIntensityImageFilter; Emmanuel Christophe committed Jan 31, 2009 71 Victor Poughon committed May 16, 2019 72 using WriterType = otb::ImageFileWriter; Emmanuel Christophe committed Jan 31, 2009 73 OTB Bot committed Apr 01, 2010 74 75 ReaderType::Pointer reader = ReaderType::New(); WriterType::Pointer writer = WriterType::New(); Victor Poughon committed Feb 20, 2019 76 HarrisFilterType::Pointer harris = HarrisFilterType::New(); OTB Bot committed Apr 01, 2010 77 RescalerType::Pointer rescaler = RescalerType::New(); Emmanuel Christophe committed Jan 31, 2009 78 OTB Bot committed Apr 01, 2010 79 80 reader->SetFileName(inputFilename); writer->SetFileName(outputFilename); Emmanuel Christophe committed Jan 31, 2009 81 OTB Bot committed Apr 01, 2010 82 harris->SetInput(reader->GetOutput()); Emmanuel Christophe committed Jan 31, 2009 83 84 85 86 87 88 89 90 91 92 93 94 // The \doxygen{otb}{HarrisImageFilter} needs some parameters to // operate. The derivative computation is performed by a // convolution with the derivative of a Gaussian kernel of // variance $\sigma_D$ (derivation scale) and // the smoothing of the image is performed by convolving with a // Gaussian kernel of variance $\sigma_I$ (integration // scale). This allows the computation of the following matrix: // // \mu(\mathbf{x},\sigma_I,\sigma_D) = \sigma_D^2 g(\sigma_I)\star // \left[\begin{array}{cc} L_x^2(\mathbf{x},\sigma_D) & // L_xL_y(\mathbf{x},\sigma_D)\\ L_xL_y(\mathbf{x},\sigma_D)& Julien Malik committed Jun 28, 2012 95 96 // L_y^2(\mathbf{x},\sigma_D) \end{array}\right] // Emmanuel Christophe committed Jan 31, 2009 97 98 // The output of the detector is $$det(\mu) - \alpha trace^2(\mu).$$ OTB Bot committed Apr 01, 2010 99 100 101 harris->SetSigmaD(SigmaD); harris->SetSigmaI(SigmaI); harris->SetAlpha(Alpha); Emmanuel Christophe committed Jan 31, 2009 102 OTB Bot committed Apr 01, 2010 103 104 rescaler->SetOutputMinimum(itk::NumericTraits::min()); rescaler->SetOutputMaximum(itk::NumericTraits::max()); Emmanuel Christophe committed Jan 31, 2009 105 OTB Bot committed Apr 01, 2010 106 107 rescaler->SetInput(harris->GetOutput()); writer->SetInput(rescaler->GetOutput()); Emmanuel Christophe committed Jan 31, 2009 108 109 writer->Update(); Jordi Inglada committed Apr 04, 2006 110 111 112 113 114 115 // Figure~\ref{fig:Harris} shows the result of applying the interest // point detector to a small patch extracted from a Spot 5 image. // \begin{figure} // \center // \includegraphics[width=0.25\textwidth]{ROISpot5.eps} // \includegraphics[width=0.25\textwidth]{ROISpot5Harris.eps} Jordi Inglada committed May 14, 2008 116 // \itkcaption[Harris Filter Application]{Result of applying the Emmanuel Christophe committed Dec 06, 2008 117 // \doxygen{otb}{HarrisImageFilter} to a Spot 5 image.} Jordi Inglada committed Apr 04, 2006 118 119 120 // \label{fig:Harris} // \end{figure} // Jordi Inglada committed Jul 03, 2006 121 // The output of the \doxygen{otb}{HarrisImageFilter} is an image OTB Bot committed Jan 31, 2011 122 // where, for each pixel, we obtain the intensity of the Jordi Inglada committed Apr 13, 2006 123 124 125 // detection. Often, the user may want to get access to the set of // points for which the output of the detector is higher than a // given threshold. This can be obtained by using the Jordi Inglada committed Jul 03, 2006 126 // \doxygen{otb}{HarrisImageToPointSetFilter}. This filter is only Jordi Inglada committed Apr 13, 2006 127 // templated over the input image type, the output being a Jordi Inglada committed Jul 03, 2006 128 // \doxygen{itk}{PointSet} with pixel type equal to the image pixel type. Jordi Inglada committed Apr 04, 2006 129 Victor Poughon committed May 16, 2019 130 using FunctionType = otb::HarrisImageToPointSetFilter; Emmanuel Christophe committed Jan 31, 2009 131 132 // We declare now the filter and a pointer to the output point set. Victor Poughon committed May 16, 2019 133 using OutputPointSetType = FunctionType::OutputPointSetType; Emmanuel Christophe committed Jan 31, 2009 134 Victor Poughon committed Feb 20, 2019 135 136 FunctionType::Pointer harrisPoints = FunctionType::New(); OutputPointSetType::Pointer pointSet = OutputPointSetType::New(); Emmanuel Christophe committed Jan 31, 2009 137 138 139 140 141 // The \doxygen{otb}{HarrisImageToPointSetFilter} takes the same // parameters as the \doxygen{otb}{HarrisImageFilter} and an // additional parameter : the threshold for the point selection. OTB Bot committed Apr 01, 2010 142 143 144 145 146 harrisPoints->SetInput(0, reader->GetOutput()); harrisPoints->SetSigmaD(SigmaD); harrisPoints->SetSigmaI(SigmaI); harrisPoints->SetAlpha(Alpha); harrisPoints->SetLowerThreshold(10); Emmanuel Christophe committed Jan 31, 2009 147 148 149 150 151 152 153 154 155 156 157 pointSet = harrisPoints->GetOutput(); harrisPoints->Update(); // We can now iterate through the obtained pointset and access // the coordinates of the points. We start by accessing the // container of the points which is encapsulated into the point // set (see section \ref{sec:PointSetSection} for more // information on using \doxygen{itk}{PointSet}s) and declaring // an iterator to it. Victor Poughon committed May 16, 2019 158 159 160 161 using ContainerType = OutputPointSetType::PointsContainer; ContainerType* pointsContainer = pointSet->GetPoints(); using IteratorType = ContainerType::Iterator; IteratorType itList = pointsContainer->Begin(); Emmanuel Christophe committed Jan 31, 2009 162 163 164 // And we get the points coordinates OTB Bot committed Apr 01, 2010 165 while (itList != pointsContainer->End()) Victor Poughon committed Feb 20, 2019 166 { Victor Poughon committed May 16, 2019 167 168 using OutputPointType = OutputPointSetType::PointType; OutputPointType pCoordinate = (itList.Value()); Guillaume Pasero committed Mar 20, 2019 169 otbLogMacro(Debug, << pCoordinate); Emmanuel Christophe committed Jan 31, 2009 170 ++itList; Victor Poughon committed Feb 20, 2019 171 } Jordi Inglada committed Apr 13, 2006 172 Jordi Inglada committed Apr 04, 2006 173 174 return EXIT_SUCCESS; }