1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
|
/*
* Copyright (C) 2005-2020 Centre National d'Etudes Spatiales (CNES)
*
* 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.
*/
/* Example usage:
./BandMathFilterExample Input/qb_RoadExtract.tif Output/RoadExtractBandMath.tif Output/qb_BandMath-pretty.jpg
*/
#include "itkMacro.h"
#include <iostream>
#include "otbImage.h"
#include "otbVectorImage.h"
#include "otbImageFileReader.h"
#include "otbImageFileWriter.h"
#include "itkUnaryFunctorImageFilter.h"
#include "itkCastImageFilter.h"
#include "otbVectorImageToImageListFilter.h"
// We start by including the required header file.
// The aim of this example is to compute the Normalized Difference Vegetation Index (NDVI)
// from a multispectral image and then apply a threshold to this
// index to extract areas containing a dense vegetation canopy.
#include "otbBandMathImageFilter.h"
int main(int argc, char* argv[])
{
if (argc != 4)
{
std::cerr << "Usage: " << argv[0] << " inputImageFile ";
std::cerr << " outputImageFile ";
std::cerr << " outputPrettyImageFile" << std::endl;
return EXIT_FAILURE;
}
// We start by the typedefs needed for reading and
// writing the images. The BandMathImageFilter class
// works with Image as input, so we need to define additional
// filters to extract each layer of the multispectral image.
using PixelType = double;
using InputImageType = otb::VectorImage<PixelType, 2>;
using OutputImageType = otb::Image<PixelType, 2>;
using ImageListType = otb::ImageList<OutputImageType>;
using VectorImageToImageListType = otb::VectorImageToImageListFilter<InputImageType, ImageListType>;
using ReaderType = otb::ImageFileReader<InputImageType>;
using WriterType = otb::ImageFileWriter<OutputImageType>;
// We can now define the type for the filter
using FilterType = otb::BandMathImageFilter<OutputImageType>;
// We instantiate the filter, the reader, and the writer
ReaderType::Pointer reader = ReaderType::New();
WriterType::Pointer writer = WriterType::New();
FilterType::Pointer filter = FilterType::New();
writer->SetInput(filter->GetOutput());
reader->SetFileName(argv[1]);
writer->SetFileName(argv[2]);
reader->UpdateOutputInformation();
// We now need to extract each band from the input VectorImage,
// it illustrates the use of the VectorImageToImageList.
// Each extracted layer is an input to the BandMathImageFilter
VectorImageToImageListType::Pointer imageList = VectorImageToImageListType::New();
imageList->SetInput(reader->GetOutput());
imageList->UpdateOutputInformation();
const unsigned int nbBands = reader->GetOutput()->GetNumberOfComponentsPerPixel();
for (unsigned int j = 0; j < nbBands; ++j)
{
filter->SetNthInput(j, imageList->GetOutput()->GetNthElement(j));
}
// Now we can define the mathematical expression to perform on the layers (b1, b2, b3, b4).
// The filter takes advantage of the parsing capabilities of the muParser library and
// allows setting the expression as on a digital calculator.
// The expression below returns 255 if the ratio (NIR-RED)/(NIR+RED) is greater than 0.4 and 0 if not.
filter->SetExpression("if((b4-b3)/(b4+b3) > 0.4, 255, 0)");
#ifdef OTB_MUPARSER_HAS_CXX_LOGICAL_OPERATORS
filter->SetExpression("((b4-b3)/(b4+b3) > 0.4) ? 255 : 0");
#else
filter->SetExpression("if((b4-b3)/(b4+b3) > 0.4, 255, 0)");
#endif
// We can now run the pipeline
writer->Update();
// The muParser library also provides the possibility to extend existing built-in functions. For example,
// you can use the OTB expression "ndvi(b3, b4)" with the filter. In this instance, the mathematical expression would be "if(ndvi(b3, b4)>0.4, 255, 0)", which
// would return the same result.
using OutputPrettyImageType = otb::Image<unsigned char, 2>;
using PrettyImageFileWriterType = otb::ImageFileWriter<OutputPrettyImageType>;
using CastImageFilterType = itk::CastImageFilter<OutputImageType, OutputPrettyImageType>;
PrettyImageFileWriterType::Pointer prettyWriter = PrettyImageFileWriterType::New();
CastImageFilterType::Pointer caster = CastImageFilterType::New();
caster->SetInput(filter->GetOutput());
prettyWriter->SetInput(caster->GetOutput());
prettyWriter->SetFileName(argv[3]);
prettyWriter->Update();
}
|
