Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: Viewshed Generator
4 : * Purpose: Viewshed Generator mainline.
5 : * Author: Tamas Szekeres <szekerest@gmail.com>
6 : *
7 : ******************************************************************************
8 : *
9 : * SPDX-License-Identifier: MIT
10 : ****************************************************************************/
11 :
12 : #include <limits>
13 :
14 : #include "commonutils.h"
15 : #include "gdal.h"
16 : #include "gdalargumentparser.h"
17 :
18 : #include "viewshed/cumulative.h"
19 : #include "viewshed/viewshed.h"
20 :
21 : namespace gdal
22 : {
23 :
24 : namespace
25 : {
26 :
27 : struct Options
28 : {
29 : viewshed::Options opts;
30 : std::string osSrcFilename;
31 : int nBandIn{1};
32 : bool bQuiet;
33 : };
34 :
35 : /// Parse arguments into options structure.
36 : ///
37 : /// \param argParser Argument parser
38 : /// \param aosArgv Command line options as a string list
39 : /// \return Command line parsed as options
40 22 : Options parseArgs(GDALArgumentParser &argParser, const CPLStringList &aosArgv)
41 : {
42 22 : Options localOpts;
43 :
44 22 : viewshed::Options &opts = localOpts.opts;
45 :
46 22 : argParser.add_output_format_argument(opts.outputFormat);
47 22 : argParser.add_argument("-ox")
48 22 : .store_into(opts.observer.x)
49 44 : .metavar("<value>")
50 22 : .help(_("The X position of the observer (in SRS units)."));
51 :
52 22 : argParser.add_argument("-oy")
53 22 : .store_into(opts.observer.y)
54 44 : .metavar("<value>")
55 22 : .help(_("The Y position of the observer (in SRS units)."));
56 :
57 22 : argParser.add_argument("-oz")
58 22 : .default_value(2)
59 22 : .store_into(opts.observer.z)
60 44 : .metavar("<value>")
61 22 : .nargs(1)
62 : .help(_("The height of the observer above the DEM surface in the "
63 22 : "height unit of the DEM."));
64 :
65 22 : argParser.add_argument("-vv")
66 22 : .default_value(255)
67 22 : .store_into(opts.visibleVal)
68 44 : .metavar("<value>")
69 22 : .nargs(1)
70 22 : .help(_("Pixel value to set for visible areas."));
71 :
72 22 : argParser.add_argument("-iv")
73 22 : .default_value(0)
74 22 : .store_into(opts.invisibleVal)
75 44 : .metavar("<value>")
76 22 : .nargs(1)
77 22 : .help(_("Pixel value to set for invisible areas."));
78 :
79 22 : argParser.add_argument("-ov")
80 22 : .default_value(0)
81 22 : .store_into(opts.outOfRangeVal)
82 44 : .metavar("<value>")
83 22 : .nargs(1)
84 : .help(
85 : _("Pixel value to set for the cells that fall outside of the range "
86 22 : "specified by the observer location and the maximum distance."));
87 :
88 22 : argParser.add_creation_options_argument(opts.creationOpts);
89 :
90 22 : argParser.add_argument("-a_nodata")
91 22 : .default_value(-1.0)
92 22 : .store_into(opts.nodataVal)
93 44 : .metavar("<value>")
94 22 : .nargs(1)
95 : .help(_("The value to be set for the cells in the output raster that "
96 22 : "have no data."));
97 :
98 22 : argParser.add_argument("-tz")
99 22 : .default_value(0.0)
100 22 : .store_into(opts.targetHeight)
101 44 : .metavar("<value>")
102 22 : .nargs(1)
103 : .help(_("The height of the target above the DEM surface in the height "
104 22 : "unit of the DEM."));
105 :
106 22 : argParser.add_argument("-md")
107 22 : .default_value(0)
108 22 : .store_into(opts.maxDistance)
109 44 : .metavar("<value>")
110 22 : .nargs(1)
111 22 : .help(_("Maximum distance from observer to compute visibility."));
112 :
113 22 : argParser.add_argument("-j")
114 22 : .default_value(3)
115 22 : .store_into(opts.numJobs)
116 44 : .metavar("<value>")
117 22 : .nargs(1)
118 : .help(_(
119 22 : "Number of relative simultaneous jobs to run in cumulative mode"));
120 :
121 : // Value for standard atmospheric refraction. See
122 : // doc/source/programs/gdal_viewshed.rst
123 22 : argParser.add_argument("-cc")
124 22 : .default_value(0.85714)
125 22 : .store_into(opts.curveCoeff)
126 44 : .metavar("<value>")
127 22 : .nargs(1)
128 : .help(_("Coefficient to consider the effect of the curvature and "
129 22 : "refraction."));
130 :
131 22 : argParser.add_argument("-b")
132 22 : .default_value(localOpts.nBandIn)
133 22 : .store_into(localOpts.nBandIn)
134 44 : .metavar("<value>")
135 22 : .nargs(1)
136 22 : .help(_("Select an input band band containing the DEM data."));
137 :
138 22 : argParser.add_argument("-om")
139 22 : .choices("NORMAL", "DEM", "GROUND", "ACCUM")
140 44 : .metavar("NORMAL|DEM|GROUND|ACCUM")
141 : .action(
142 38 : [&into = opts.outputMode](const std::string &value)
143 : {
144 8 : if (EQUAL(value.c_str(), "DEM"))
145 1 : into = viewshed::OutputMode::DEM;
146 7 : else if (EQUAL(value.c_str(), "GROUND"))
147 1 : into = viewshed::OutputMode::Ground;
148 6 : else if (EQUAL(value.c_str(), "ACCUM"))
149 1 : into = viewshed::OutputMode::Cumulative;
150 : else
151 5 : into = viewshed::OutputMode::Normal;
152 22 : })
153 22 : .nargs(1)
154 22 : .help(_("Sets what information the output contains."));
155 :
156 22 : argParser.add_argument("-os")
157 22 : .default_value(10)
158 22 : .store_into(opts.observerSpacing)
159 44 : .metavar("<value>")
160 22 : .nargs(1)
161 22 : .help(_("Spacing between observer cells when using cumulative mode."));
162 :
163 22 : argParser.add_quiet_argument(&localOpts.bQuiet);
164 :
165 22 : argParser.add_argument("src_filename")
166 22 : .store_into(localOpts.osSrcFilename)
167 22 : .metavar("<src_filename>");
168 :
169 22 : argParser.add_argument("dst_filename")
170 22 : .store_into(opts.outputFilename)
171 22 : .metavar("<dst_filename>");
172 :
173 : try
174 : {
175 22 : argParser.parse_args(aosArgv);
176 : }
177 2 : catch (const std::exception &err)
178 : {
179 2 : argParser.display_error_and_usage(err);
180 2 : std::exit(1);
181 : }
182 19 : return localOpts;
183 : }
184 :
185 : /// Validate specified options.
186 : ///
187 : /// \param localOpts Options to validate
188 : /// \param argParser Argument parser
189 19 : void validateArgs(Options &localOpts, const GDALArgumentParser &argParser)
190 : {
191 19 : viewshed::Options &opts = localOpts.opts;
192 :
193 19 : if (opts.maxDistance < 0)
194 : {
195 0 : CPLError(CE_Failure, CPLE_AppDefined,
196 : "Max distance must be non-negative.");
197 0 : exit(2);
198 : }
199 :
200 19 : if (opts.outputFormat.empty())
201 : {
202 : opts.outputFormat =
203 12 : GetOutputDriverForRaster(opts.outputFilename.c_str());
204 12 : if (opts.outputFormat.empty())
205 : {
206 0 : exit(2);
207 : }
208 : }
209 :
210 19 : if (opts.outputMode != viewshed::OutputMode::Cumulative)
211 : {
212 54 : for (const char *opt : {"-os", "-j"})
213 36 : if (argParser.is_used(opt))
214 : {
215 0 : std::string err = "Option " + std::string(opt) +
216 0 : " can only be used in cumulative mode.";
217 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", err.c_str());
218 0 : exit(2);
219 : }
220 : }
221 :
222 19 : if (opts.outputMode == viewshed::OutputMode::Cumulative)
223 : {
224 6 : for (const char *opt : {"-ox", "-oy", "-vv", "-iv", "-md"})
225 5 : if (argParser.is_used(opt))
226 : {
227 0 : std::string err = "Option " + std::string(opt) +
228 0 : " can't be used in cumulative mode.";
229 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", err.c_str());
230 0 : exit(2);
231 : }
232 : }
233 : else
234 : {
235 51 : for (const char *opt : {"-ox", "-oy"})
236 35 : if (!argParser.is_used(opt))
237 : {
238 : std::string err =
239 4 : "Option " + std::string(opt) + " is required.";
240 2 : CPLError(CE_Failure, CPLE_AppDefined, "%s", err.c_str());
241 2 : exit(2);
242 : }
243 : }
244 :
245 : // For double values that are out of range for byte raster output,
246 : // set to zero. Values less than zero are sentinel as NULL nodata.
247 31 : if (opts.outputMode == viewshed::OutputMode::Normal &&
248 14 : opts.nodataVal > std::numeric_limits<uint8_t>::max())
249 0 : opts.nodataVal = 0;
250 17 : }
251 :
252 : // Adjust the coefficient of curvature for non-earth SRS.
253 : /// \param curveCoeff Current curve coefficient
254 : /// \param hSrcDS Source dataset
255 : /// \return Adjusted curve coefficient.
256 12 : double adjustCurveCoeff(double curveCoeff, GDALDatasetH hSrcDS)
257 : {
258 : const OGRSpatialReference *poSRS =
259 12 : GDALDataset::FromHandle(hSrcDS)->GetSpatialRef();
260 12 : if (poSRS)
261 : {
262 : OGRErr eSRSerr;
263 10 : const double dfSemiMajor = poSRS->GetSemiMajor(&eSRSerr);
264 10 : if (eSRSerr != OGRERR_FAILURE &&
265 10 : fabs(dfSemiMajor - SRS_WGS84_SEMIMAJOR) >
266 : 0.05 * SRS_WGS84_SEMIMAJOR)
267 : {
268 1 : curveCoeff = 1.0;
269 1 : CPLDebug("gdal_viewshed",
270 : "Using -cc=1.0 as a non-Earth CRS has been detected");
271 : }
272 : }
273 12 : return curveCoeff;
274 : }
275 :
276 : } // unnamed namespace
277 : } // namespace gdal
278 :
279 : /************************************************************************/
280 : /* main() */
281 : /************************************************************************/
282 :
283 22 : MAIN_START(argc, argv)
284 : {
285 : using namespace gdal;
286 :
287 22 : EarlySetConfigOptions(argc, argv);
288 :
289 22 : GDALAllRegister();
290 :
291 22 : argc = GDALGeneralCmdLineProcessor(argc, &argv, 0);
292 37 : CPLStringList aosArgv;
293 22 : aosArgv.Assign(argv, /* bTakeOwnership= */ true);
294 22 : if (argc < 1)
295 0 : std::exit(-argc);
296 :
297 59 : GDALArgumentParser argParser(aosArgv[0], /* bForBinary=*/true);
298 :
299 : argParser.add_description(
300 22 : _("Calculates a viewshed raster from an input raster DEM."));
301 :
302 : argParser.add_epilog(_("For more details, consult "
303 22 : "https://gdal.org/programs/gdal_viewshed.html"));
304 :
305 22 : Options localOpts = parseArgs(argParser, aosArgv);
306 19 : viewshed::Options &opts = localOpts.opts;
307 :
308 19 : validateArgs(localOpts, argParser);
309 :
310 : /* -------------------------------------------------------------------- */
311 : /* Open source raster file. */
312 : /* -------------------------------------------------------------------- */
313 : GDALDatasetH hSrcDS =
314 17 : GDALOpen(localOpts.osSrcFilename.c_str(), GA_ReadOnly);
315 17 : if (hSrcDS == nullptr)
316 1 : exit(2);
317 :
318 16 : GDALRasterBandH hBand = GDALGetRasterBand(hSrcDS, localOpts.nBandIn);
319 16 : if (hBand == nullptr)
320 : {
321 1 : CPLError(CE_Failure, CPLE_AppDefined,
322 : "Band %d does not exist on dataset.", localOpts.nBandIn);
323 1 : exit(2);
324 : }
325 :
326 15 : if (!argParser.is_used("-cc"))
327 12 : opts.curveCoeff = adjustCurveCoeff(opts.curveCoeff, hSrcDS);
328 :
329 : /* -------------------------------------------------------------------- */
330 : /* Invoke. */
331 : /* -------------------------------------------------------------------- */
332 :
333 : GDALDatasetH hDstDS;
334 :
335 : bool bSuccess;
336 15 : if (opts.outputMode == viewshed::OutputMode::Cumulative)
337 : {
338 1 : viewshed::Cumulative oViewshed(opts);
339 1 : bSuccess = oViewshed.run(localOpts.osSrcFilename,
340 1 : localOpts.bQuiet ? GDALDummyProgress
341 : : GDALTermProgress);
342 : }
343 : else
344 : {
345 28 : viewshed::Viewshed oViewshed(opts);
346 14 : bSuccess = oViewshed.run(hBand, localOpts.bQuiet ? GDALDummyProgress
347 : : GDALTermProgress);
348 14 : hDstDS = GDALDataset::FromHandle(oViewshed.output().release());
349 14 : GDALClose(hSrcDS);
350 14 : if (GDALClose(hDstDS) != CE_None)
351 0 : bSuccess = false;
352 : }
353 :
354 15 : GDALDestroyDriverManager();
355 15 : OGRCleanupAll();
356 :
357 15 : return bSuccess ? 0 : 1;
358 : }
359 :
360 0 : MAIN_END
|