otbEndmemberNumberEstimation.cxx 6.5 KB
Newer Older
1
/*
2
 * Copyright (C) 2005-2020 Centre National d'Etudes Spatiales (CNES)
3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
 *
 * 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.
 */

#include "otbWrapperApplication.h"
#include "otbWrapperApplicationFactory.h"

24
#include "otbStreamingStatisticsVectorImageFilter.h"
25 26
#include "otbEigenvalueLikelihoodMaximisation.h"
#include "otbVirtualDimensionality.h"
27

28 29 30 31 32 33 34 35 36
namespace otb
{
namespace Wrapper
{

class EndmemberNumberEstimation : public Application
{
public:
  /** Standard class typedefs. */
37 38 39 40
  typedef EndmemberNumberEstimation     Self;
  typedef Application                   Superclass;
  typedef itk::SmartPointer<Self>       Pointer;
  typedef itk::SmartPointer<const Self> ConstPointer;
41 42 43 44 45 46

  /** Standard macro */
  itkNewMacro(Self);

  itkTypeMacro(EndmemberNumberEstimation, otb::Application);

47
  typedef otb::StreamingStatisticsVectorImageFilter<FloatVectorImageType, double> StreamingStatisticsVectorImageFilterType;
48 49
  typedef otb::VirtualDimensionality<double>            VirtualDimensionalityType;
  typedef otb::EigenvalueLikelihoodMaximisation<double> EigenvalueLikelihoodMaximisationType;
50

51 52 53 54
private:
  void DoInit() override
  {
    SetName("EndmemberNumberEstimation");
55
    SetDescription("Estimate the number of endmembers in a hyperspectral image");
56 57

    // Documentation
58 59 60 61 62
    SetDocLongDescription(
        "Estimate the number of endmembers "
        "in a hyperspectral image. First, compute statistics on the image and then "
        "apply an endmember number estimation algorithm using these statistics. Two "
        "algorithms are available:\n\n"
63

64 65
        "1. Virtual Dimensionality (HFC-VD) [1][2]\n"
        "2. Eigenvalue Likelihood Maximization (ELM) [3][4]\n\n"
66

67
        "The application then returns the estimated number of endmembers.\n\n"
68

69 70 71
        "[1] C.-I. Chang and Q. Du, Estimation of number of spectrally distinct signal "
        "sources in hyperspectral imagery, IEEE Transactions on Geoscience and Remote "
        "Sensing, vol. 43, no. 3, mar 2004.\n\n"
72

73 74 75 76
        "[2] J. Wang and C.-I. Chang, Applications of independent component analysis "
        "in endmember extraction and abundance quantification for hyperspectral imagery"
        ", IEEE Transactions on Geoscience and Remote Sensing, vol. 44, no. 9, pp. "
        "2601-1616, sep 2006.\n\n"
77

78 79
        "[3] Unsupervised Endmember Extraction of Martian Hyperspectral Images, B.Luo, "
        "J. Chanussot, S. Dout\'e and X. Ceamanos, IEEE Whispers 2009, Grenoble France, 2009\n\n"
80

81 82 83
        "[4] Unsupervised classification of hyperspectral images by using "
        "linear unmixing algorithm Luo, B. and Chanussot, J., IEEE Int. Conf. On Image"
        "Processing(ICIP) 2009, Cairo, Egypte, 2009");
84

85 86 87 88 89 90
    SetDocLimitations("None");
    SetDocAuthors("OTB-Team");
    SetDocSeeAlso("VertexComponentAnalysis, HyperspectralUnmixing");

    AddDocTag(Tags::Hyperspectral);

91 92
    AddParameter(ParameterType_InputImage, "in", "Input Image Filename");
    SetParameterDescription("in", "The hyperspectral data cube input");
93

94
    AddParameter(ParameterType_Int, "number", "Number of endmembers");
95 96 97
    SetParameterDescription("number", "The output estimated number of endmembers");
    SetParameterRole("number", Role_Output);

98 99
    AddParameter(ParameterType_Choice, "algo", "Unmixing algorithm");
    SetParameterDescription("algo", "The algorithm to use for the estimation");
100
    AddChoice("algo.elm", "Eigenvalue Likelihood Maximization");
101
    SetParameterDescription("algo.elm", "Eigenvalue Likelihood Maximization algorithm");
102
    AddChoice("algo.vd", "Virtual Dimensionality");
103
    SetParameterDescription("algo.vd", "HFC Virtual Dimensionality algorithm");
104

105
    AddParameter(ParameterType_Float, "algo.vd.far", "False alarm rate");
Cédric Traizet's avatar
Cédric Traizet committed
106 107
    SetMinimumParameterFloatValue("algo.vd.far", 0);
    SetMaximumParameterFloatValue("algo.vd.far", 1);
108 109
    SetDefaultParameterFloat("algo.vd.far", 1.0E-3);
    SetParameterDescription("algo.vd.far", "False alarm rate for the virtual dimensionality algorithm");
110

111
    AddRAMParameter();
112

113 114
    // Doc example parameter settings
    SetDocExampleParameterValue("in", "cupriteSubHsi.tif");
115 116
    SetDocExampleParameterValue("algo", "vd");
    SetDocExampleParameterValue("algo.vd.far", "1.0E-3");
117 118 119 120 121 122 123 124 125 126 127

    SetOfficialDocLink();
  }

  void DoUpdateParameters() override
  {
    // Nothing to do here : all parameters are independent
  }

  void DoExecute() override
  {
Cédric Traizet's avatar
Cédric Traizet committed
128 129 130
    // Load input image
    auto inputImage = GetParameterImage("in");

131
    otbAppLogINFO("Computing statistics on input image...");
132
    auto statisticsFilter = StreamingStatisticsVectorImageFilterType::New();
Cédric Traizet's avatar
Cédric Traizet committed
133
    statisticsFilter->SetInput(inputImage);
134
    AddProcess(statisticsFilter->GetStreamer(), "Statistic estimation step");
135

136
    statisticsFilter->Update();
Cédric Traizet's avatar
Cédric Traizet committed
137 138

    auto correlationMatrix = statisticsFilter->GetCorrelation().GetVnlMatrix();
139 140
    auto covarianceMatrix  = statisticsFilter->GetCovariance().GetVnlMatrix();
    auto numberOfPixels    = inputImage->GetLargestPossibleRegion().GetNumberOfPixels();
Cédric Traizet's avatar
Cédric Traizet committed
141

142 143 144 145
    int               numberOfEndmembers = 0;
    const std::string algorithm          = GetParameterString("algo");
    if (algorithm == "elm")
    {
146
      otbAppLogINFO("Estimation algorithm: Eigenvalue Likelihood Maximization.");
147
      auto elm = EigenvalueLikelihoodMaximisationType::New();
Cédric Traizet's avatar
Cédric Traizet committed
148 149 150
      elm->SetCovariance(covarianceMatrix);
      elm->SetCorrelation(correlationMatrix);
      elm->SetNumberOfPixels(numberOfPixels);
151 152
      elm->Compute();
      numberOfEndmembers = elm->GetNumberOfEndmembers();
153 154 155
    }
    else if (algorithm == "vd")
    {
156
      otbAppLogINFO("Estimation algorithm: HFC Virtual Dimensionality.");
157
      auto vd = VirtualDimensionalityType::New();
Cédric Traizet's avatar
Cédric Traizet committed
158 159 160
      vd->SetCovariance(covarianceMatrix);
      vd->SetCorrelation(correlationMatrix);
      vd->SetNumberOfPixels(numberOfPixels);
161
      vd->SetFAR(GetParameterFloat("algo.vd.far"));
162 163
      vd->Compute();
      numberOfEndmembers = vd->GetNumberOfEndmembers();
164
    }
165
    SetParameterInt("number", numberOfEndmembers);
166 167 168 169 170 171
  }
};
}
}

OTB_APPLICATION_EXPORT(otb::Wrapper::EndmemberNumberEstimation)