OGS
createIntermediateRasters.cpp File Reference
#include <tclap/CmdLine.h>
#include <memory>
#include <string>
#include "BaseLib/FileTools.h"
#include "BaseLib/Logging.h"
#include "BaseLib/MPI.h"
#include "BaseLib/TCLAPArguments.h"
#include "GeoLib/IO/AsciiRasterInterface.h"
#include "GeoLib/Raster.h"
#include "InfoLib/GitInfo.h"
Include dependency graph for createIntermediateRasters.cpp:

Go to the source code of this file.

Functions

int main (int argc, char *argv[])

Function Documentation

◆ main()

int main ( int argc,
char * argv[] )

Definition at line 17 of file createIntermediateRasters.cpp.

18{
19 TCLAP::CmdLine cmd(
20 "Takes two DEMs located at the exact same spatial position (but at "
21 "different elevation) and calculates n raster DEMs located at "
22 "equidistant intervals between them (i.e. for n=1, one new raster "
23 "located precisely in the middle will be created).\n\n"
24 "OpenGeoSys-6 software, version " +
26 ".\n"
27 "Copyright (c) 2012-2026, OpenGeoSys Community "
28 "(http://www.opengeosys.org)",
30 TCLAP::ValueArg<std::size_t> number_arg(
31 "n", "number",
32 "number of rasters to be calculated, "
33 "(min = 0)",
34 false, 1, "int");
35 cmd.add(number_arg);
36 TCLAP::ValueArg<std::string> output_arg("o", "output-file",
37 "Output (.asc). Raster output file",
38 true, "", "OUTPUT_FILE");
39 cmd.add(output_arg);
40 TCLAP::ValueArg<std::string> input2_arg(
41 "", "file2", "Input (.asc). Second DEM-raster input file", true, "",
42 "INPUT_FILE");
43 cmd.add(input2_arg);
44 TCLAP::ValueArg<std::string> input1_arg(
45 "", "file1", "Input (.asc) First DEM-raster input file", true, "",
46 "INPUT_FILE");
47 cmd.add(input1_arg);
48 auto log_level_arg = BaseLib::makeLogLevelArg();
49 cmd.add(log_level_arg);
50
51 cmd.parse(argc, argv);
52
53 BaseLib::MPI::Setup mpi_setup(argc, argv);
54 BaseLib::initOGSLogger(log_level_arg.getValue());
55
56 std::unique_ptr<GeoLib::Raster> dem1(
57 FileIO::AsciiRasterInterface::readRaster(input1_arg.getValue()));
58 std::unique_ptr<GeoLib::Raster> dem2(
59 FileIO::AsciiRasterInterface::readRaster(input2_arg.getValue()));
60
61 if (dem1 == nullptr || dem2 == nullptr)
62 {
63 return 1;
64 }
65
66 GeoLib::RasterHeader const& h1 = dem1->getHeader();
67 GeoLib::RasterHeader const& h2 = dem2->getHeader();
68
69 bool errors_found(false);
70 if (std::abs(h1.origin[0] - h2.origin[0]) > (h1.cell_size / 100.0))
71 {
72 ERR("Origin x-coordinate is not the same in both raster files.\n");
73 errors_found = true;
74 }
75 if (std::abs(h1.origin[1] - h2.origin[1]) > (h1.cell_size / 100.0))
76 {
77 ERR("Origin y-coordinate is not the same in both raster files.\n");
78 errors_found = true;
79 }
80 if (h1.cell_size != h2.cell_size)
81 {
82 ERR("Cellsize is not the same in both raster files.\n");
83 errors_found = true;
84 }
85 if (h1.n_cols != h2.n_cols)
86 {
87 ERR("Raster width is not the same in both raster files.\n");
88 errors_found = true;
89 }
90 if (h1.n_rows != h2.n_rows)
91 {
92 ERR("Raster height is not the same in both raster files.\n");
93 errors_found = true;
94 }
95
96 if (errors_found)
97 {
98 return 2;
99 }
100
101 std::size_t const n = number_arg.getValue();
102 std::vector<std::vector<double>> raster;
103 for (std::size_t i = 0; i < n; ++i)
104 {
105 std::vector<double> r;
106 r.reserve(h1.n_cols * h1.n_rows);
107 raster.push_back(r);
108 }
109
110 auto it2 = dem2->begin();
111 for (auto it1 = dem1->begin(); it1 != dem1->end(); ++it1)
112 {
113 if (it2 == dem2->end())
114 {
115 ERR("Error: File 2 is shorter than File 1.");
116 return 3;
117 }
118 if (*it1 == h1.no_data || *it2 == h2.no_data)
119 {
120 for (std::size_t i = 0; i < n; ++i)
121 {
122 raster[i].push_back(h1.no_data);
123 }
124 }
125 else
126 {
127 double const min = std::min(*it1, *it2);
128 double const max = std::max(*it1, *it2);
129 double const step = (max - min) / static_cast<double>(n + 1);
130 for (std::size_t i = 0; i < n; ++i)
131 {
132 raster[i].push_back(max - ((i + 1) * step));
133 }
134 }
135 it2++;
136 }
137 if (it2 != dem2->end())
138 {
139 ERR("Error: File 1 is shorter than File 2.");
140 return 3;
141 }
142
143 std::string const filename = output_arg.getValue();
144 for (std::size_t i = 0; i < n; ++i)
145 {
146 std::string const basename = BaseLib::dropFileExtension(filename);
147 std::string const ext = BaseLib::getFileExtension(filename);
148
149 GeoLib::Raster r(h1, raster[i].begin(), raster[i].end());
151 r, basename + std::to_string(i) + ext);
152 INFO("Layer {:d} written.", i + 1);
153 }
154 return EXIT_SUCCESS;
155}
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:28
void ERR(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:40
static void writeRasterAsASC(GeoLib::Raster const &raster, std::string const &file_name)
Writes an Esri asc-file.
static GeoLib::Raster * readRaster(std::string const &fname)
Class Raster is used for managing raster data.
Definition Raster.h:39
TCLAP::ValueArg< std::string > makeLogLevelArg()
void initOGSLogger(std::string const &log_level)
Definition Logging.cpp:56
std::string getFileExtension(const std::string &path)
std::string dropFileExtension(std::string const &filename)
GITINFOLIB_EXPORT const std::string ogs_version
Contains the relevant information when storing a geoscientific raster data.
Definition Raster.h:18
MathLib::Point3d origin
Definition Raster.h:22
std::size_t n_cols
Definition Raster.h:19
std::size_t n_rows
Definition Raster.h:20

References GeoLib::RasterHeader::cell_size, BaseLib::dropFileExtension(), ERR(), BaseLib::getFileExtension(), INFO(), BaseLib::initOGSLogger(), BaseLib::makeLogLevelArg(), GeoLib::RasterHeader::n_cols, GeoLib::RasterHeader::n_rows, GeoLib::RasterHeader::no_data, GitInfoLib::GitInfo::ogs_version, GeoLib::RasterHeader::origin, FileIO::AsciiRasterInterface::readRaster(), and FileIO::AsciiRasterInterface::writeRasterAsASC().