Line data    Source code

       1             : /******************************************************************************
       2             :  *
       3             :  * Project:  GDAL
       4             :  * Purpose:  gdal "raster contour" subcommand
       5             :  * Author:   Alessandro Pasotti <elpaso at itopen dot it>
       6             :  *
       7             :  ******************************************************************************
       8             :  * Copyright (c) 2025, Alessandro Pasotti <elpaso at itopen dot it>
       9             :  *
      10             :  * SPDX-License-Identifier: MIT
      11             :  ****************************************************************************/
      12             : 
      13             : #include <cmath>
      14             : 
      15             : #include "gdalalg_raster_contour.h"
      16             : 
      17             : #include "cpl_conv.h"
      18             : #include "gdal_priv.h"
      19             : #include "gdal_utils.h"
      20             : #include "gdal_alg.h"
      21             : #include "gdal_utils_priv.h"
      22             : 
      23             : //! @cond Doxygen_Suppress
      24             : 
      25             : #ifndef _
      26             : #define _(x) (x)
      27             : #endif
      28             : 
      29             : /************************************************************************/
      30             : /*          GDALRasterContourAlgorithm::GDALRasterContourAlgorithm()    */
      31             : /************************************************************************/
      32             : 
      33          19 : GDALRasterContourAlgorithm::GDALRasterContourAlgorithm()
      34             :     : GDALAlgorithm(NAME, DESCRIPTION, HELP_URL), m_outputLayerName("contour"),
      35          19 :       m_elevAttributeName(""), m_amin(""), m_amax(""), m_levels{}
      36             : {
      37             : 
      38          19 :     AddProgressArg();
      39          19 :     AddOutputFormatArg(&m_outputFormat)
      40             :         .AddMetadataItem(GAAMDI_REQUIRED_CAPABILITIES,
      41          57 :                          {GDAL_DCAP_VECTOR, GDAL_DCAP_CREATE});
      42          19 :     AddOpenOptionsArg(&m_openOptions);
      43          19 :     AddInputFormatsArg(&m_inputFormats)
      44          38 :         .AddMetadataItem(GAAMDI_REQUIRED_CAPABILITIES, {GDAL_DCAP_RASTER});
      45          19 :     AddInputDatasetArg(&m_inputDataset, GDAL_OF_RASTER);
      46          19 :     AddOutputDatasetArg(&m_outputDataset, GDAL_OF_VECTOR);
      47          19 :     AddCreationOptionsArg(&m_creationOptions);
      48          19 :     AddLayerCreationOptionsArg(&m_layerCreationOptions);
      49             : 
      50             :     // gdal_contour specific options
      51          19 :     AddBandArg(&m_band).SetDefault(1);
      52          19 :     AddLayerNameArg(&m_outputLayerName).AddAlias("nln");
      53             :     AddArg("elevation-name", 0, _("Name of the elevation field"),
      54          19 :            &m_elevAttributeName);
      55          19 :     AddArg("min-name", 0, _("Name of the minimum elevation field"), &m_amin);
      56          19 :     AddArg("max-name", 0, _("Name of the maximum elevation field"), &m_amax);
      57          19 :     AddArg("3d", 0, _("Force production of 3D vectors instead of 2D"), &m_3d);
      58             : 
      59             :     AddArg("src-nodata", 0, _("Input pixel value to treat as 'nodata'"),
      60          19 :            &m_sNodata);
      61          38 :     AddArg("interval", 0, _("Elevation interval between contours"), &m_interval)
      62          38 :         .SetMutualExclusionGroup("levels")
      63          19 :         .SetMinValueExcluded(0);
      64          38 :     AddArg("levels", 0, _("List of contour levels"), &m_levels)
      65          19 :         .SetMutualExclusionGroup("levels");
      66             :     AddArg("exp-base", 'e', _("Base for exponential contour level generation"),
      67          38 :            &m_expBase)
      68          19 :         .SetMutualExclusionGroup("levels");
      69          38 :     AddArg("offset", 0, _("Offset to apply to contour levels"), &m_offset)
      70          19 :         .AddAlias("off");
      71             :     AddArg("polygonize", 'p', _("Create polygons instead of lines"),
      72          19 :            &m_polygonize);
      73             :     AddArg("group-transactions", 0,
      74             :            _("Group n features per transaction (default 100 000)"),
      75          38 :            &m_groupTransactions)
      76          19 :         .SetMinValueIncluded(0);
      77          19 :     AddOverwriteArg(&m_overwrite);
      78          19 : }
      79             : 
      80             : /************************************************************************/
      81             : /*                  GDALRasterContourAlgorithm::RunImpl()               */
      82             : /************************************************************************/
      83             : 
      84          13 : bool GDALRasterContourAlgorithm::RunImpl(GDALProgressFunc pfnProgress,
      85             :                                          void *pProgressData)
      86             : {
      87             : 
      88          13 :     CPLErrorReset();
      89             : 
      90          13 :     CPLAssert(m_inputDataset.GetDatasetRef());
      91          13 :     CPLAssert(!m_outputDataset.GetDatasetRef());
      92             : 
      93          26 :     CPLStringList aosOptions;
      94          13 :     if (!m_outputFormat.empty())
      95             :     {
      96           0 :         aosOptions.AddString("-of");
      97           0 :         aosOptions.AddString(m_outputFormat);
      98             :     }
      99             : 
     100          14 :     for (const auto &co : m_creationOptions)
     101             :     {
     102           1 :         aosOptions.AddString("-co");
     103           1 :         aosOptions.AddString(co);
     104             :     }
     105             : 
     106          14 :     for (const auto &co : m_layerCreationOptions)
     107             :     {
     108           1 :         aosOptions.AddString("-lco");
     109           1 :         aosOptions.AddString(co);
     110             :     }
     111             : 
     112          13 :     if (m_band > 0)
     113             :     {
     114          13 :         aosOptions.AddString("-b");
     115          13 :         aosOptions.AddString(CPLSPrintf("%d", m_band));
     116             :     }
     117          13 :     if (!m_elevAttributeName.empty())
     118             :     {
     119           5 :         aosOptions.AddString("-a");
     120           5 :         aosOptions.AddString(m_elevAttributeName);
     121             :     }
     122          13 :     if (!m_amin.empty())
     123             :     {
     124           6 :         aosOptions.AddString("-amin");
     125           6 :         aosOptions.AddString(m_amin);
     126             :     }
     127          13 :     if (!m_amax.empty())
     128             :     {
     129           6 :         aosOptions.AddString("-amax");
     130           6 :         aosOptions.AddString(m_amax);
     131             :     }
     132          13 :     if (m_3d)
     133             :     {
     134           0 :         aosOptions.AddString("-3d");
     135             :     }
     136          13 :     if (!std::isnan(m_sNodata))
     137             :     {
     138           1 :         aosOptions.AddString("-snodata");
     139           1 :         aosOptions.AddString(CPLSPrintf("%.16g", m_sNodata));
     140             :     }
     141          13 :     if (m_levels.size() > 0)
     142             :     {
     143          12 :         for (const auto &level : m_levels)
     144             :         {
     145           8 :             aosOptions.AddString("-fl");
     146           8 :             aosOptions.AddString(level);
     147             :         }
     148             :     }
     149          13 :     if (!std::isnan(m_interval))
     150             :     {
     151           7 :         aosOptions.AddString("-i");
     152           7 :         aosOptions.AddString(CPLSPrintf("%.16g", m_interval));
     153             :     }
     154          13 :     if (m_expBase > 0)
     155             :     {
     156           1 :         aosOptions.AddString("-e");
     157           1 :         aosOptions.AddString(CPLSPrintf("%d", m_expBase));
     158             :     }
     159          13 :     if (!std::isnan(m_offset))
     160             :     {
     161           1 :         aosOptions.AddString("-off");
     162           1 :         aosOptions.AddString(CPLSPrintf("%.16g", m_offset));
     163             :     }
     164          13 :     if (m_polygonize)
     165             :     {
     166           5 :         aosOptions.AddString("-p");
     167             :     }
     168          13 :     if (!m_outputLayerName.empty())
     169             :     {
     170          13 :         aosOptions.AddString("-nln");
     171          13 :         aosOptions.AddString(m_outputLayerName);
     172             :     }
     173             : 
     174             :     // Check that one of --interval, --levels, --exp-base is specified
     175          13 :     if (m_levels.size() == 0 && std::isnan(m_interval) && m_expBase == 0)
     176             :     {
     177           1 :         ReportError(
     178             :             CE_Failure, CPLE_AppDefined,
     179             :             "One of 'interval', 'levels', 'exp-base' must be specified.");
     180           1 :         return false;
     181             :     }
     182             : 
     183             :     // Check that interval is not negative
     184          12 :     if (!std::isnan(m_interval) && m_interval < 0)
     185             :     {
     186           0 :         ReportError(CE_Failure, CPLE_AppDefined,
     187             :                     "Interval must be a positive number.");
     188           0 :         return false;
     189             :     }
     190             : 
     191          12 :     aosOptions.AddString(m_inputDataset.GetName());
     192          12 :     aosOptions.AddString(m_outputDataset.GetName());
     193             : 
     194          24 :     GDALContourOptionsForBinary optionsForBinary;
     195             :     std::unique_ptr<GDALContourOptions, decltype(&GDALContourOptionsFree)>
     196             :         psOptions{GDALContourOptionsNew(aosOptions.List(), &optionsForBinary),
     197          24 :                   GDALContourOptionsFree};
     198             : 
     199          12 :     if (!psOptions)
     200             :     {
     201           0 :         return false;
     202             :     }
     203             : 
     204          12 :     GDALDatasetH hSrcDS{m_inputDataset.GetDatasetRef()};
     205          12 :     GDALRasterBandH hBand{nullptr};
     206          12 :     GDALDatasetH hDstDS{m_outputDataset.GetDatasetRef()};
     207          12 :     OGRLayerH hLayer{nullptr};
     208          12 :     char **papszStringOptions = nullptr;
     209             : 
     210             :     CPLErr eErr =
     211          12 :         GDALContourProcessOptions(psOptions.get(), &papszStringOptions, &hSrcDS,
     212             :                                   &hBand, &hDstDS, &hLayer);
     213             : 
     214          12 :     if (eErr == CE_None)
     215             :     {
     216          12 :         eErr = GDALContourGenerateEx(hBand, hLayer, papszStringOptions,
     217             :                                      pfnProgress, pProgressData);
     218             :     }
     219             : 
     220          12 :     CSLDestroy(papszStringOptions);
     221             : 
     222          12 :     auto poDstDS = GDALDataset::FromHandle(hDstDS);
     223          12 :     m_outputDataset.Set(std::unique_ptr<GDALDataset>(poDstDS));
     224             : 
     225          12 :     return eErr == CE_None;
     226             : }
     227             : 
     228             : //! @endcond