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          27 : GDALRasterContourAlgorithm::GDALRasterContourAlgorithm(bool standaloneStep)
      34             :     : GDALPipelineStepAlgorithm(
      35             :           NAME, DESCRIPTION, HELP_URL,
      36           0 :           ConstructorOptions()
      37          27 :               .SetStandaloneStep(standaloneStep)
      38          54 :               .SetOutputFormatCreateCapability(GDAL_DCAP_CREATE))
      39             : {
      40          27 :     m_outputLayerName = "contour";
      41             : 
      42          27 :     AddProgressArg();
      43          27 :     if (standaloneStep)
      44             :     {
      45          20 :         AddOutputFormatArg(&m_format).AddMetadataItem(
      46          60 :             GAAMDI_REQUIRED_CAPABILITIES, {GDAL_DCAP_VECTOR, GDAL_DCAP_CREATE});
      47          20 :         AddOpenOptionsArg(&m_openOptions);
      48          20 :         AddInputFormatsArg(&m_inputFormats)
      49          40 :             .AddMetadataItem(GAAMDI_REQUIRED_CAPABILITIES, {GDAL_DCAP_RASTER});
      50          20 :         AddInputDatasetArg(&m_inputDataset, GDAL_OF_RASTER);
      51          20 :         AddOutputDatasetArg(&m_outputDataset, GDAL_OF_VECTOR);
      52          20 :         AddCreationOptionsArg(&m_creationOptions);
      53          20 :         AddLayerCreationOptionsArg(&m_layerCreationOptions);
      54          20 :         AddLayerNameArg(&m_outputLayerName).AddAlias("nln");
      55          20 :         AddOverwriteArg(&m_overwrite);
      56             :     }
      57             : 
      58             :     // gdal_contour specific options
      59          27 :     AddBandArg(&m_band).SetDefault(1);
      60             : 
      61             :     AddArg("elevation-name", 0, _("Name of the elevation field"),
      62          27 :            &m_elevAttributeName);
      63          27 :     AddArg("min-name", 0, _("Name of the minimum elevation field"), &m_amin);
      64          27 :     AddArg("max-name", 0, _("Name of the maximum elevation field"), &m_amax);
      65          27 :     AddArg("3d", 0, _("Force production of 3D vectors instead of 2D"), &m_3d);
      66             : 
      67             :     AddArg("src-nodata", 0, _("Input pixel value to treat as 'nodata'"),
      68          27 :            &m_sNodata);
      69          54 :     AddArg("interval", 0, _("Elevation interval between contours"), &m_interval)
      70          54 :         .SetMutualExclusionGroup("levels")
      71          27 :         .SetMinValueExcluded(0);
      72          54 :     AddArg("levels", 0, _("List of contour levels"), &m_levels)
      73          27 :         .SetMutualExclusionGroup("levels");
      74             :     AddArg("exp-base", 'e', _("Base for exponential contour level generation"),
      75          54 :            &m_expBase)
      76          27 :         .SetMutualExclusionGroup("levels");
      77          54 :     AddArg("offset", 0, _("Offset to apply to contour levels"), &m_offset)
      78          27 :         .AddAlias("off");
      79             :     AddArg("polygonize", 'p', _("Create polygons instead of lines"),
      80          27 :            &m_polygonize);
      81             :     AddArg("group-transactions", 0,
      82             :            _("Group n features per transaction (default 100 000)"),
      83          54 :            &m_groupTransactions)
      84          27 :         .SetMinValueIncluded(0);
      85          27 : }
      86             : 
      87             : /************************************************************************/
      88             : /*                  GDALRasterContourAlgorithm::RunImpl()               */
      89             : /************************************************************************/
      90             : 
      91          13 : bool GDALRasterContourAlgorithm::RunImpl(GDALProgressFunc pfnProgress,
      92             :                                          void *pProgressData)
      93             : {
      94          13 :     GDALPipelineStepRunContext stepCtxt;
      95          13 :     stepCtxt.m_pfnProgress = pfnProgress;
      96          13 :     stepCtxt.m_pProgressData = pProgressData;
      97          13 :     return RunPreStepPipelineValidations() && RunStep(stepCtxt);
      98             : }
      99             : 
     100             : /************************************************************************/
     101             : /*                  GDALRasterContourAlgorithm::RunStep()               */
     102             : /************************************************************************/
     103             : 
     104          14 : bool GDALRasterContourAlgorithm::RunStep(GDALPipelineStepRunContext &ctxt)
     105             : {
     106          14 :     CPLErrorReset();
     107             : 
     108          14 :     auto poSrcDS = m_inputDataset[0].GetDatasetRef();
     109          14 :     CPLAssert(poSrcDS);
     110             : 
     111          14 :     CPLAssert(!m_outputDataset.GetDatasetRef());
     112             : 
     113          28 :     CPLStringList aosOptions;
     114             : 
     115          28 :     std::string outputFilename;
     116          14 :     if (m_standaloneStep)
     117             :     {
     118          13 :         outputFilename = m_outputDataset.GetName();
     119          13 :         if (!m_format.empty())
     120             :         {
     121           0 :             aosOptions.AddString("-of");
     122           0 :             aosOptions.AddString(m_format.c_str());
     123             :         }
     124             : 
     125          14 :         for (const auto &co : m_creationOptions)
     126             :         {
     127           1 :             aosOptions.push_back("-co");
     128           1 :             aosOptions.push_back(co.c_str());
     129             :         }
     130             : 
     131          14 :         for (const auto &co : m_layerCreationOptions)
     132             :         {
     133           1 :             aosOptions.push_back("-lco");
     134           1 :             aosOptions.push_back(co.c_str());
     135             :         }
     136             :     }
     137             :     else
     138             :     {
     139           1 :         if (GetGDALDriverManager()->GetDriverByName("GPKG"))
     140             :         {
     141           1 :             aosOptions.AddString("-of");
     142           1 :             aosOptions.AddString("GPKG");
     143             : 
     144           1 :             outputFilename = CPLGenerateTempFilenameSafe("_contour.gpkg");
     145             :         }
     146             :         else
     147             :         {
     148           0 :             aosOptions.AddString("-of");
     149           0 :             aosOptions.AddString("MEM");
     150             :         }
     151             :     }
     152             : 
     153          14 :     if (m_band > 0)
     154             :     {
     155          14 :         aosOptions.AddString("-b");
     156          14 :         aosOptions.AddString(CPLSPrintf("%d", m_band));
     157             :     }
     158          14 :     if (!m_elevAttributeName.empty())
     159             :     {
     160           5 :         aosOptions.AddString("-a");
     161           5 :         aosOptions.AddString(m_elevAttributeName);
     162             :     }
     163          14 :     if (!m_amin.empty())
     164             :     {
     165           6 :         aosOptions.AddString("-amin");
     166           6 :         aosOptions.AddString(m_amin);
     167             :     }
     168          14 :     if (!m_amax.empty())
     169             :     {
     170           6 :         aosOptions.AddString("-amax");
     171           6 :         aosOptions.AddString(m_amax);
     172             :     }
     173          14 :     if (m_3d)
     174             :     {
     175           0 :         aosOptions.AddString("-3d");
     176             :     }
     177          14 :     if (!std::isnan(m_sNodata))
     178             :     {
     179           1 :         aosOptions.AddString("-snodata");
     180           1 :         aosOptions.AddString(CPLSPrintf("%.16g", m_sNodata));
     181             :     }
     182          14 :     if (m_levels.size() > 0)
     183             :     {
     184          12 :         for (const auto &level : m_levels)
     185             :         {
     186           8 :             aosOptions.AddString("-fl");
     187           8 :             aosOptions.AddString(level);
     188             :         }
     189             :     }
     190          14 :     if (!std::isnan(m_interval))
     191             :     {
     192           8 :         aosOptions.AddString("-i");
     193           8 :         aosOptions.AddString(CPLSPrintf("%.16g", m_interval));
     194             :     }
     195          14 :     if (m_expBase > 0)
     196             :     {
     197           1 :         aosOptions.AddString("-e");
     198           1 :         aosOptions.AddString(CPLSPrintf("%d", m_expBase));
     199             :     }
     200          14 :     if (!std::isnan(m_offset))
     201             :     {
     202           1 :         aosOptions.AddString("-off");
     203           1 :         aosOptions.AddString(CPLSPrintf("%.16g", m_offset));
     204             :     }
     205          14 :     if (m_polygonize)
     206             :     {
     207           5 :         aosOptions.AddString("-p");
     208             :     }
     209          14 :     if (!m_outputLayerName.empty())
     210             :     {
     211          14 :         aosOptions.AddString("-nln");
     212          14 :         aosOptions.AddString(m_outputLayerName);
     213             :     }
     214             : 
     215             :     // Check that one of --interval, --levels, --exp-base is specified
     216          14 :     if (m_levels.size() == 0 && std::isnan(m_interval) && m_expBase == 0)
     217             :     {
     218           1 :         ReportError(
     219             :             CE_Failure, CPLE_AppDefined,
     220             :             "One of 'interval', 'levels', 'exp-base' must be specified.");
     221           1 :         return false;
     222             :     }
     223             : 
     224             :     // Check that interval is not negative
     225          13 :     if (!std::isnan(m_interval) && m_interval < 0)
     226             :     {
     227           0 :         ReportError(CE_Failure, CPLE_AppDefined,
     228             :                     "Interval must be a positive number.");
     229           0 :         return false;
     230             :     }
     231             : 
     232          13 :     aosOptions.AddString(m_inputDataset[0].GetName());
     233          13 :     aosOptions.AddString(outputFilename);
     234             : 
     235          13 :     bool bRet = false;
     236          26 :     GDALContourOptionsForBinary optionsForBinary;
     237             :     std::unique_ptr<GDALContourOptions, decltype(&GDALContourOptionsFree)>
     238             :         psOptions{GDALContourOptionsNew(aosOptions.List(), &optionsForBinary),
     239          13 :                   GDALContourOptionsFree};
     240          13 :     if (psOptions)
     241             :     {
     242          13 :         GDALDatasetH hSrcDS{poSrcDS};
     243          13 :         GDALRasterBandH hBand{nullptr};
     244          13 :         GDALDatasetH hDstDS{m_outputDataset.GetDatasetRef()};
     245          13 :         OGRLayerH hLayer{nullptr};
     246          13 :         char **papszStringOptions = nullptr;
     247             : 
     248          13 :         bRet = GDALContourProcessOptions(psOptions.get(), &papszStringOptions,
     249             :                                          &hSrcDS, &hBand, &hDstDS,
     250             :                                          &hLayer) == CE_None;
     251             : 
     252          13 :         if (bRet)
     253             :         {
     254          13 :             bRet = GDALContourGenerateEx(hBand, hLayer, papszStringOptions,
     255             :                                          ctxt.m_pfnProgress,
     256             :                                          ctxt.m_pProgressData) == CE_None;
     257             :         }
     258             : 
     259          13 :         CSLDestroy(papszStringOptions);
     260             : 
     261          13 :         auto poDstDS = GDALDataset::FromHandle(hDstDS);
     262          13 :         if (bRet)
     263             :         {
     264          13 :             bRet = poDstDS != nullptr;
     265             :         }
     266          13 :         if (poDstDS && !m_standaloneStep && !outputFilename.empty())
     267             :         {
     268           1 :             poDstDS->MarkSuppressOnClose();
     269           1 :             if (bRet)
     270           1 :                 bRet = poDstDS->FlushCache() == CE_None;
     271           1 :             VSIUnlink(outputFilename.c_str());
     272             :         }
     273          13 :         m_outputDataset.Set(std::unique_ptr<GDALDataset>(poDstDS));
     274             :     }
     275             : 
     276          13 :     return bRet;
     277             : }
     278             : 
     279             : GDALRasterContourAlgorithmStandalone::~GDALRasterContourAlgorithmStandalone() =
     280             :     default;
     281             : 
     282             : //! @endcond