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