OTB  9.0.0
Orfeo Toolbox
BasicFilters/HillShadingExample.cxx
/*
* Copyright (C) 2005-2022 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:
./HillShadingExample Output/HillShadingExample.png Output/HillShadingColorExample.png 6.5 45.5 500 500 0.002 -0.002 Input/DEM/srtm_directory
*/
#include "itkScalarToRGBColormapImageFilter.h"
#include "itkMultiplyImageFilter.h"
#include "itkShiftScaleImageFilter.h"
#include "otbWorldFile.h"
int main(int argc, char* argv[])
{
if (argc < 10)
{
std::cout << argv[0] << " <output_filename> <output_color_filename> "
<< " <Longitude Output Origin point> <Latitude Output Origin point>"
<< " <X Output Size> <Y Output size>"
<< " <X Spacing> <Y Spacing> <DEM folder path>"
<< " [Projection Ref]" << std::endl;
return EXIT_FAILURE;
}
using PixelType = double;
using UCharPixelType = unsigned char;
using RGBPixelType = itk::RGBPixel<UCharPixelType>;
using ImageType = otb::Image<PixelType, 2>;
using RGBImageType = otb::Image<RGBPixelType, 2>;
using ScalarImageType = otb::Image<UCharPixelType, 2>;
using ScalarWriterType = otb::ImageFileWriter<ScalarImageType>;
ScalarWriterType::Pointer writer = ScalarWriterType::New();
writer->SetFileName(argv[1]);
WriterType::Pointer writer2 = WriterType::New();
writer2->SetFileName(argv[2]);
using DEMToImageGeneratorType = otb::DEMToImageGenerator<ImageType>;
DEMToImageGeneratorType::Pointer demToImage = DEMToImageGeneratorType::New();
using SizeType = DEMToImageGeneratorType::SizeType;
using SpacingType = DEMToImageGeneratorType::SpacingType;
using PointType = DEMToImageGeneratorType::PointType;
PointType origin;
origin[0] = ::atof(argv[3]);
origin[1] = ::atof(argv[4]);
demToImage->SetOutputOrigin(origin);
SizeType size;
size[0] = ::atoi(argv[5]);
size[1] = ::atoi(argv[6]);
demToImage->SetOutputSize(size);
SpacingType spacing;
spacing[0] = ::atof(argv[7]);
spacing[1] = ::atof(argv[8]);
demToImage->SetOutputSpacing(spacing);
double res = 0;
if (argc > 10)
{
demToImage->SetOutputProjectionRef(argv[10]);
res = spacing[0];
}
else
{
// Compute the resolution (Vincenty formula)
double lon1 = origin[0];
double lon2 = origin[0] + size[0] * spacing[0];
double lat1 = origin[1];
double lat2 = origin[1] + size[1] * spacing[1];
double R = 6371; // km
double d = std::acos(std::sin(lat1) * std::sin(lat2) + std::cos(lat1) * std::cos(lat2) * std::cos(lon2 - lon1)) * R;
res = d / std::sqrt(2.0);
}
// After generating the DEM image as in the DEMToImageGenerator example, you can declare
// the hill shading mechanism. The hill shading is implemented as a functor doing some
// operations in its neighborhood. A convenient filter called \doxygen{otb}{HillShadingFilter}
// is defined around this mechanism.
using HillShadingFilterType = otb::HillShadingFilter<ImageType, ImageType>;
HillShadingFilterType::Pointer hillShading = HillShadingFilterType::New();
hillShading->SetRadius(1);
hillShading->SetInput(demToImage->GetOutput());
hillShading->GetFunctor().SetXRes(res);
hillShading->GetFunctor().SetYRes(res);
using RescalerType = itk::ShiftScaleImageFilter<ImageType, ScalarImageType>;
RescalerType::Pointer rescaler = RescalerType::New();
rescaler->SetScale(255.0);
rescaler->SetInput(hillShading->GetOutput());
writer->SetInput(rescaler->GetOutput());
using ColorMapFilterType = itk::ScalarToRGBColormapImageFilter<ImageType, RGBImageType>;
ColorMapFilterType::Pointer colormapper = ColorMapFilterType::New();
colormapper->UseInputImageExtremaForScalingOff();
ColorMapFunctorType::Pointer colormap = ColorMapFunctorType::New();
colormap->SetMinimumInputValue(0);
colormap->SetMaximumInputValue(4000);
colormapper->SetColormap(colormap);
colormapper->SetInput(demToImage->GetOutput());
using MultiplyFilterType =
itk::BinaryFunctorImageFilter<RGBImageType, ImageType, RGBImageType, otb::Functor::HillShadeModulationFunctor<RGBPixelType, PixelType, RGBPixelType>>;
MultiplyFilterType::Pointer multiply = MultiplyFilterType::New();
multiply->SetInput1(colormapper->GetOutput());
multiply->SetInput2(hillShading->GetOutput());
writer2->SetInput(multiply->GetOutput());
writer->Update();
writer2->Update();
worldFile->SetLonOrigin(origin[0]);
worldFile->SetLatOrigin(origin[1]);
worldFile->SetLonSpacing(spacing[0]);
worldFile->SetLatSpacing(spacing[1]);
worldFile->SetImageFilename(argv[1]);
worldFile->Update();
worldFile->SetImageFilename(argv[2]);
worldFile->Update();
}
otb::DEMHandler::GetInstance
static DEMHandler & GetInstance()
otb::WorldFile::Pointer
itk::SmartPointer< Self > Pointer
Definition: otbWorldFile.h:50
otbDEMToImageGenerator.h
otb::Image
Creation of an "otb" image which contains metadata.
Definition: otbImage.h:89
otb::HillShadingFilter
Filter to compute the lambertian of a surface given as an image.
Definition: otbHillShadingFilter.h:42
otb::DEMToImageGenerator
Class to generate an image from DEM data.
Definition: otbDEMToImageGenerator.h:51
otbWorldFile.h
main
int main(int ac, char *av[])
Definition: otbTestMain.h:88
otb::Functor::ReliefColormapFunctor
Function object which maps a scalar value into an RGB colormap value for relief representation.
Definition: otbReliefColormapFunctor.h:46
otb::ImageFileWriter
Writes image data to a single file with streaming process.
Definition: otbImageFileWriter.h:66
otbImageFileWriter.h
otbReliefColormapFunctor.h
otbImageFileReader.h
otb::WorldFile::New
static Pointer New()
otbHillShadingFilter.h
otb::DEMHandler::OpenDEMDirectory
void OpenDEMDirectory(std::string DEMDirectory)