Coverage Report

Created: 2026-02-14 06:52

next uncovered line (L), next uncovered region (R), next uncovered branch (B)
/src/gdal/apps/gdalalg_raster_viewshed.cpp
Line
Count
Source
1
/******************************************************************************
2
 *
3
 * Project:  GDAL
4
 * Purpose:  gdal "raster viewshed" subcommand
5
 * Author:   Even Rouault <even dot rouault at spatialys.com>
6
 *
7
 ******************************************************************************
8
 * Copyright (c) 2025, Even Rouault <even dot rouault at spatialys.com>
9
 *
10
 * SPDX-License-Identifier: MIT
11
 ****************************************************************************/
12
13
#include "gdalalg_raster_viewshed.h"
14
15
#include "cpl_conv.h"
16
#include "cpl_vsi_virtual.h"
17
18
#include "commonutils.h"
19
#include "gdal_priv.h"
20
21
#include "viewshed/cumulative.h"
22
#include "viewshed/viewshed.h"
23
24
#include <algorithm>
25
26
//! @cond Doxygen_Suppress
27
28
#ifndef _
29
0
#define _(x) (x)
30
#endif
31
32
/************************************************************************/
33
/*      GDALRasterViewshedAlgorithm::GDALRasterViewshedAlgorithm()      */
34
/************************************************************************/
35
36
GDALRasterViewshedAlgorithm::GDALRasterViewshedAlgorithm(bool standaloneStep)
37
0
    : GDALRasterPipelineNonNativelyStreamingAlgorithm(NAME, DESCRIPTION,
38
0
                                                      HELP_URL, standaloneStep)
39
0
{
40
0
    AddArg("position", 'p', _("Observer position"), &m_observerPos)
41
0
        .AddAlias("pos")
42
0
        .SetMetaVar("<X,Y> or <X,Y,H>")
43
0
        .SetMinCount(2)
44
0
        .SetMaxCount(3)
45
0
        .SetRepeatedArgAllowed(false);
46
0
    AddArg("height", 'z', _("Observer height"), &m_opts.observer.z);
47
48
0
    auto &sdFilenameArg =
49
0
        AddArg("sd-filename", 0, _("Filename of standard-deviation raster"),
50
0
               &m_sdFilename, GDAL_OF_RASTER);
51
0
    SetAutoCompleteFunctionForFilename(sdFilenameArg, GDAL_OF_RASTER);
52
53
0
    AddArg("target-height", 0,
54
0
           _("Height of the target above the DEM surface in the height unit of "
55
0
             "the DEM."),
56
0
           &m_opts.targetHeight)
57
0
        .SetDefault(m_opts.targetHeight);
58
0
    AddArg("mode", 0, _("Sets what information the output contains."),
59
0
           &m_outputMode)
60
0
        .SetChoices("normal", "DEM", "ground", "cumulative")
61
0
        .SetDefault(m_outputMode);
62
63
0
    AddArg("max-distance", 0,
64
0
           _("Maximum distance from observer to compute visibility. It is also "
65
0
             "used to clamp the extent of the output raster."),
66
0
           &m_opts.maxDistance)
67
0
        .SetMinValueIncluded(0);
68
0
    AddArg("min-distance", 0,
69
0
           _("Mask all cells less than this distance from the observer. Must "
70
0
             "be less "
71
0
             "than 'max-distance'."),
72
0
           &m_opts.minDistance)
73
0
        .SetMinValueIncluded(0);
74
75
0
    AddArg("start-angle", 0,
76
0
           _("Mask all cells outside of the arc ('start-angle', 'end-angle'). "
77
0
             "Clockwise degrees "
78
0
             "from north. Also used to clamp the extent of the output raster."),
79
0
           &m_opts.startAngle)
80
0
        .SetMinValueIncluded(0)
81
0
        .SetMaxValueExcluded(360);
82
0
    AddArg("end-angle", 0,
83
0
           _("Mask all cells outside of the arc ('start-angle', 'end-angle'). "
84
0
             "Clockwise degrees "
85
0
             "from north. Also used to clamp the extent of the output raster."),
86
0
           &m_opts.endAngle)
87
0
        .SetMinValueIncluded(0)
88
0
        .SetMaxValueExcluded(360);
89
90
0
    AddArg("high-pitch", 0,
91
0
           _("Mark all cells out-of-range where the observable height would be "
92
0
             "higher than the "
93
0
             "'high-pitch' angle from the observer. Degrees from horizontal - "
94
0
             "positive is up. "
95
0
             "Must be greater than 'low-pitch'."),
96
0
           &m_opts.highPitch)
97
0
        .SetMaxValueIncluded(90)
98
0
        .SetMinValueExcluded(-90);
99
0
    AddArg("low-pitch", 0,
100
0
           _("Bound observable height to be no lower than the 'low-pitch' "
101
0
             "angle from the observer. "
102
0
             "Degrees from horizontal - positive is up. Must be less than "
103
0
             "'high-pitch'."),
104
0
           &m_opts.lowPitch)
105
0
        .SetMaxValueExcluded(90)
106
0
        .SetMinValueIncluded(-90);
107
108
0
    AddArg("curvature-coefficient", 0,
109
0
           _("Coefficient to consider the effect of the curvature and "
110
0
             "refraction."),
111
0
           &m_opts.curveCoeff)
112
0
        .SetMinValueIncluded(0);
113
114
0
    AddBandArg(&m_band).SetDefault(m_band);
115
0
    AddArg("visible-value", 0, _("Pixel value to set for visible areas"),
116
0
           &m_opts.visibleVal)
117
0
        .SetDefault(m_opts.visibleVal)
118
0
        .SetMinValueIncluded(0)
119
0
        .SetMaxValueIncluded(255);
120
0
    AddArg("invisible-value", 0, _("Pixel value to set for invisible areas"),
121
0
           &m_opts.invisibleVal)
122
0
        .SetDefault(m_opts.invisibleVal)
123
0
        .SetMinValueIncluded(0)
124
0
        .SetMaxValueIncluded(255);
125
0
    AddArg("maybe-visible-value", 0,
126
0
           _("Pixel value to set for potentially visible areas"),
127
0
           &m_opts.maybeVisibleVal)
128
0
        .SetDefault(m_opts.maybeVisibleVal)
129
0
        .SetMinValueIncluded(0)
130
0
        .SetMaxValueIncluded(255);
131
0
    AddArg("out-of-range-value", 0,
132
0
           _("Pixel value to set for the cells that fall outside of the range "
133
0
             "specified by the observer location and the maximum distance"),
134
0
           &m_opts.outOfRangeVal)
135
0
        .SetDefault(m_opts.outOfRangeVal)
136
0
        .SetMinValueIncluded(0)
137
0
        .SetMaxValueIncluded(255);
138
0
    AddArg("dst-nodata", 0,
139
0
           _("The value to be set for the cells in the output raster that have "
140
0
             "no data."),
141
0
           &m_opts.nodataVal)
142
0
        .SetMinValueIncluded(0)
143
0
        .SetMaxValueIncluded(255);
144
0
    AddArg("observer-spacing", 0, _("Cell Spacing between observers"),
145
0
           &m_opts.observerSpacing)
146
0
        .SetDefault(m_opts.observerSpacing)
147
0
        .SetMinValueIncluded(1);
148
149
0
    m_numThreadsStr = std::to_string(m_numThreads);
150
0
    AddNumThreadsArg(&m_numThreads, &m_numThreadsStr);
151
0
}
152
153
/************************************************************************/
154
/*                GDALRasterViewshedAlgorithm::RunStep()                */
155
/************************************************************************/
156
157
bool GDALRasterViewshedAlgorithm::RunStep(GDALPipelineStepRunContext &ctxt)
158
0
{
159
0
    auto pfnProgress = ctxt.m_pfnProgress;
160
0
    auto pProgressData = ctxt.m_pProgressData;
161
0
    auto poSrcDS = m_inputDataset[0].GetDatasetRef();
162
0
    CPLAssert(poSrcDS);
163
0
    CPLAssert(!m_outputDataset.GetDatasetRef());
164
165
0
    GDALRasterBandH sdBand = nullptr;
166
0
    if (auto sdDataset = m_sdFilename.GetDatasetRef())
167
0
    {
168
0
        if (sdDataset->GetRasterCount() == 0)
169
0
        {
170
0
            ReportError(
171
0
                CE_Failure, CPLE_AppDefined,
172
0
                "The standard deviation dataset must have one raster band");
173
0
            return false;
174
0
        }
175
0
        sdBand = GDALRasterBand::FromHandle(sdDataset->GetRasterBand(1));
176
0
    }
177
178
0
    if (GetArg("height")->IsExplicitlySet())
179
0
    {
180
0
        if (m_observerPos.size() == 3)
181
0
        {
182
0
            ReportError(CE_Failure, CPLE_AppDefined,
183
0
                        "Height can't be specified in both 'position' and "
184
0
                        "'height' arguments");
185
0
            return false;
186
0
        }
187
0
    }
188
189
0
    if (m_observerPos.size())
190
0
    {
191
0
        m_opts.observer.x = m_observerPos[0];
192
0
        m_opts.observer.y = m_observerPos[1];
193
0
        if (m_observerPos.size() == 3)
194
0
            m_opts.observer.z = m_observerPos[2];
195
0
        else
196
0
            m_opts.observer.z = 2;
197
0
    }
198
199
0
    if (!GetArg("curvature-coefficient")->IsExplicitlySet())
200
0
    {
201
0
        m_opts.curveCoeff = gdal::viewshed::adjustCurveCoeff(
202
0
            m_opts.curveCoeff, GDALDataset::ToHandle(poSrcDS));
203
0
    }
204
205
0
    if (m_outputMode == "normal")
206
0
        m_opts.outputMode = gdal::viewshed::OutputMode::Normal;
207
0
    else if (m_outputMode == "DEM")
208
0
        m_opts.outputMode = gdal::viewshed::OutputMode::DEM;
209
0
    else if (m_outputMode == "ground")
210
0
        m_opts.outputMode = gdal::viewshed::OutputMode::Ground;
211
0
    else if (m_outputMode == "cumulative")
212
0
        m_opts.outputMode = gdal::viewshed::OutputMode::Cumulative;
213
214
0
    m_opts.numJobs = static_cast<uint8_t>(std::clamp(m_numThreads, 0, 255));
215
216
0
    m_opts.outputFilename =
217
0
        CPLGenerateTempFilenameSafe(
218
0
            CPLGetBasenameSafe(poSrcDS->GetDescription()).c_str()) +
219
0
        ".tif";
220
0
    m_opts.outputFormat = "GTiff";
221
222
0
    if (m_opts.outputMode == gdal::viewshed::OutputMode::Cumulative)
223
0
    {
224
0
        static const std::vector<std::string> badArgs{
225
0
            "visible-value", "invisible-value", "max-distance",
226
0
            "min-distance",  "start-angle",     "end-angle",
227
0
            "low-pitch",     "high-pitch",      "position"};
228
229
0
        for (const auto &arg : badArgs)
230
0
            if (GetArg(arg)->IsExplicitlySet())
231
0
            {
232
0
                std::string err =
233
0
                    "Option '" + arg + "' can't be used in cumulative mode.";
234
0
                ReportError(CE_Failure, CPLE_AppDefined, "%s", err.c_str());
235
0
                return false;
236
0
            }
237
238
0
        auto poSrcDriver = poSrcDS->GetDriver();
239
0
        if (EQUAL(poSrcDS->GetDescription(), "") || !poSrcDriver ||
240
0
            EQUAL(poSrcDriver->GetDescription(), "MEM"))
241
0
        {
242
0
            ReportError(
243
0
                CE_Failure, CPLE_AppDefined,
244
0
                "In cumulative mode, the input dataset must be opened by name");
245
0
            return false;
246
0
        }
247
0
        gdal::viewshed::Cumulative oViewshed(m_opts);
248
0
        const bool bSuccess = oViewshed.run(
249
0
            m_inputDataset[0].GetName().c_str(),
250
0
            pfnProgress ? pfnProgress : GDALDummyProgress, pProgressData);
251
0
        if (bSuccess)
252
0
        {
253
0
            m_outputDataset.Set(std::unique_ptr<GDALDataset>(
254
0
                GDALDataset::Open(m_opts.outputFilename.c_str(),
255
0
                                  GDAL_OF_RASTER | GDAL_OF_VERBOSE_ERROR,
256
0
                                  nullptr, nullptr, nullptr)));
257
0
        }
258
0
    }
259
0
    else
260
0
    {
261
0
        static const std::vector<std::string> badArgs{
262
0
            "observer-spacing", GDAL_ARG_NAME_NUM_THREADS};
263
0
        for (const auto &arg : badArgs)
264
0
            if (GetArg(arg)->IsExplicitlySet())
265
0
            {
266
0
                std::string err =
267
0
                    "Option '" + arg + "' can't be used in standard mode.";
268
0
                ReportError(CE_Failure, CPLE_AppDefined, "%s", err.c_str());
269
0
                return false;
270
0
            }
271
0
        static const std::vector<std::string> goodArgs{"position"};
272
0
        for (const auto &arg : goodArgs)
273
0
            if (!GetArg(arg)->IsExplicitlySet())
274
0
            {
275
0
                std::string err =
276
0
                    "Option '" + arg + "' must be specified in standard mode.";
277
0
                ReportError(CE_Failure, CPLE_AppDefined, "%s", err.c_str());
278
0
                return false;
279
0
            }
280
281
0
        gdal::viewshed::Viewshed oViewshed(m_opts);
282
0
        const bool bSuccess = oViewshed.run(
283
0
            GDALRasterBand::ToHandle(poSrcDS->GetRasterBand(m_band)), sdBand,
284
0
            pfnProgress ? pfnProgress : GDALDummyProgress, pProgressData);
285
0
        if (bSuccess)
286
0
        {
287
0
            m_outputDataset.Set(oViewshed.output());
288
0
        }
289
0
    }
290
291
0
    auto poOutDS = m_outputDataset.GetDatasetRef();
292
0
    if (poOutDS && poOutDS->GetDescription()[0])
293
0
    {
294
        // In file systems that allow it (all but Windows...), we want to
295
        // delete the temporary file as soon as soon as possible after
296
        // having open it, so that if someone kills the process there are
297
        // no temp files left over. If that unlink() doesn't succeed
298
        // (on Windows), then the file will eventually be deleted when
299
        // poTmpDS is cleaned due to MarkSuppressOnClose().
300
0
        VSIUnlink(poOutDS->GetDescription());
301
0
        poOutDS->MarkSuppressOnClose();
302
0
    }
303
304
0
    return poOutDS != nullptr;
305
0
}
306
307
GDALRasterViewshedAlgorithmStandalone::
308
0
    ~GDALRasterViewshedAlgorithmStandalone() = default;
309
310
//! @endcond