LCOV - code coverage report
Current view: top level - frmts/raw - lcpdataset.cpp (source / functions) Hit Total Coverage
Test: gdal_filtered.info Lines: 702 824 85.2 %
Date: 2025-10-01 17:07:58 Functions: 12 12 100.0 %

          Line data    Source code
       1             : /******************************************************************************
       2             :  *
       3             :  * Project:  LCP Driver
       4             :  * Purpose:  FARSITE v.4 Landscape file (.lcp) reader for GDAL
       5             :  * Author:   Chris Toney
       6             :  *
       7             :  ******************************************************************************
       8             :  * Copyright (c) 2008, Chris Toney
       9             :  * Copyright (c) 2008-2011, Even Rouault <even dot rouault at spatialys.com>
      10             :  * Copyright (c) 2013, Kyle Shannon <kyle at pobox dot com>
      11             :  *
      12             :  * SPDX-License-Identifier: MIT
      13             :  ****************************************************************************/
      14             : 
      15             : #include "cpl_port.h"
      16             : #include "cpl_string.h"
      17             : #include "gdal_frmts.h"
      18             : #include "gdal_priv.h"
      19             : #include "ogr_spatialref.h"
      20             : #include "rawdataset.h"
      21             : 
      22             : #include <algorithm>
      23             : #include <limits>
      24             : 
      25             : constexpr size_t LCP_HEADER_SIZE = 7316;
      26             : constexpr int LCP_MAX_BANDS = 10;
      27             : constexpr int LCP_MAX_PATH = 256;
      28             : constexpr int LCP_MAX_DESC = 512;
      29             : constexpr int LCP_MAX_CLASSES = 100;
      30             : 
      31             : /************************************************************************/
      32             : /* ==================================================================== */
      33             : /*                              LCPDataset                              */
      34             : /* ==================================================================== */
      35             : /************************************************************************/
      36             : 
      37             : class LCPDataset final : public RawDataset
      38             : {
      39             :     VSILFILE *fpImage;  // image data file.
      40             :     char pachHeader[LCP_HEADER_SIZE];
      41             : 
      42             :     CPLString osPrjFilename{};
      43             :     OGRSpatialReference m_oSRS{};
      44             : 
      45             :     static CPLErr ClassifyBandData(GDALRasterBand *poBand, GInt32 &nNumClasses,
      46             :                                    GInt32 *panClasses);
      47             : 
      48             :     CPL_DISALLOW_COPY_ASSIGN(LCPDataset)
      49             : 
      50             :     CPLErr Close() override;
      51             : 
      52             :   public:
      53             :     LCPDataset();
      54             :     ~LCPDataset() override;
      55             : 
      56             :     char **GetFileList(void) override;
      57             : 
      58             :     CPLErr GetGeoTransform(GDALGeoTransform &gt) const override;
      59             : 
      60             :     static int Identify(GDALOpenInfo *);
      61             :     static GDALDataset *Open(GDALOpenInfo *);
      62             :     static GDALDataset *CreateCopy(const char *pszFilename,
      63             :                                    GDALDataset *poSrcDS, int bStrict,
      64             :                                    char **papszOptions,
      65             :                                    GDALProgressFunc pfnProgress,
      66             :                                    void *pProgressData);
      67             : 
      68           3 :     const OGRSpatialReference *GetSpatialRef() const override
      69             :     {
      70           3 :         return &m_oSRS;
      71             :     }
      72             : };
      73             : 
      74             : /************************************************************************/
      75             : /*                            LCPDataset()                              */
      76             : /************************************************************************/
      77             : 
      78          75 : LCPDataset::LCPDataset() : fpImage(nullptr)
      79             : {
      80          75 :     memset(pachHeader, 0, sizeof(pachHeader));
      81          75 : }
      82             : 
      83             : /************************************************************************/
      84             : /*                            ~LCPDataset()                             */
      85             : /************************************************************************/
      86             : 
      87         150 : LCPDataset::~LCPDataset()
      88             : 
      89             : {
      90          75 :     LCPDataset::Close();
      91         150 : }
      92             : 
      93             : /************************************************************************/
      94             : /*                              Close()                                 */
      95             : /************************************************************************/
      96             : 
      97         132 : CPLErr LCPDataset::Close()
      98             : {
      99         132 :     CPLErr eErr = CE_None;
     100         132 :     if (nOpenFlags != OPEN_FLAGS_CLOSED)
     101             :     {
     102          75 :         if (LCPDataset::FlushCache(true) != CE_None)
     103           0 :             eErr = CE_Failure;
     104             : 
     105          75 :         if (fpImage)
     106             :         {
     107          75 :             if (VSIFCloseL(fpImage) != 0)
     108             :             {
     109           0 :                 CPLError(CE_Failure, CPLE_FileIO, "I/O error");
     110           0 :                 eErr = CE_Failure;
     111             :             }
     112             :         }
     113             : 
     114          75 :         if (GDALPamDataset::Close() != CE_None)
     115           0 :             eErr = CE_Failure;
     116             :     }
     117         132 :     return eErr;
     118             : }
     119             : 
     120             : /************************************************************************/
     121             : /*                          GetGeoTransform()                           */
     122             : /************************************************************************/
     123             : 
     124           2 : CPLErr LCPDataset::GetGeoTransform(GDALGeoTransform &gt) const
     125             : {
     126           2 :     double dfEast = 0.0;
     127           2 :     double dfWest = 0.0;
     128           2 :     double dfNorth = 0.0;
     129           2 :     double dfSouth = 0.0;
     130           2 :     double dfCellX = 0.0;
     131           2 :     double dfCellY = 0.0;
     132             : 
     133           2 :     memcpy(&dfEast, pachHeader + 4172, sizeof(double));
     134           2 :     memcpy(&dfWest, pachHeader + 4180, sizeof(double));
     135           2 :     memcpy(&dfNorth, pachHeader + 4188, sizeof(double));
     136           2 :     memcpy(&dfSouth, pachHeader + 4196, sizeof(double));
     137           2 :     memcpy(&dfCellX, pachHeader + 4208, sizeof(double));
     138           2 :     memcpy(&dfCellY, pachHeader + 4216, sizeof(double));
     139           2 :     CPL_LSBPTR64(&dfEast);
     140           2 :     CPL_LSBPTR64(&dfWest);
     141           2 :     CPL_LSBPTR64(&dfNorth);
     142           2 :     CPL_LSBPTR64(&dfSouth);
     143           2 :     CPL_LSBPTR64(&dfCellX);
     144           2 :     CPL_LSBPTR64(&dfCellY);
     145             : 
     146           2 :     gt[0] = dfWest;
     147           2 :     gt[3] = dfNorth;
     148           2 :     gt[1] = dfCellX;
     149           2 :     gt[2] = 0.0;
     150             : 
     151           2 :     gt[4] = 0.0;
     152           2 :     gt[5] = -1 * dfCellY;
     153             : 
     154           2 :     return CE_None;
     155             : }
     156             : 
     157             : /************************************************************************/
     158             : /*                              Identify()                              */
     159             : /************************************************************************/
     160             : 
     161       58762 : int LCPDataset::Identify(GDALOpenInfo *poOpenInfo)
     162             : 
     163             : {
     164             :     /* -------------------------------------------------------------------- */
     165             :     /*      Verify that this is a FARSITE v.4 LCP file                      */
     166             :     /* -------------------------------------------------------------------- */
     167       58762 :     if (poOpenInfo->nHeaderBytes < 50)
     168       54937 :         return FALSE;
     169             : 
     170             :     /* check if first three fields have valid data */
     171        3825 :     if ((CPL_LSBSINT32PTR(poOpenInfo->pabyHeader) != 20 &&
     172        3813 :          CPL_LSBSINT32PTR(poOpenInfo->pabyHeader) != 21) ||
     173         166 :         (CPL_LSBSINT32PTR(poOpenInfo->pabyHeader + 4) != 20 &&
     174         142 :          CPL_LSBSINT32PTR(poOpenInfo->pabyHeader + 4) != 21) ||
     175         166 :         (CPL_LSBSINT32PTR(poOpenInfo->pabyHeader + 8) < -90 ||
     176         168 :          CPL_LSBSINT32PTR(poOpenInfo->pabyHeader + 8) > 90))
     177             :     {
     178        3657 :         return FALSE;
     179             :     }
     180             : 
     181             : /* -------------------------------------------------------------------- */
     182             : /*      Check file extension                                            */
     183             : /* -------------------------------------------------------------------- */
     184             : #ifndef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION
     185         168 :     const char *pszFileExtension = poOpenInfo->osExtension.c_str();
     186         168 :     if (!EQUAL(pszFileExtension, "lcp"))
     187             :     {
     188           0 :         return FALSE;
     189             :     }
     190             : #endif
     191             : 
     192         168 :     return TRUE;
     193             : }
     194             : 
     195             : /************************************************************************/
     196             : /*                            GetFileList()                             */
     197             : /************************************************************************/
     198             : 
     199          38 : char **LCPDataset::GetFileList()
     200             : 
     201             : {
     202          38 :     char **papszFileList = GDALPamDataset::GetFileList();
     203             : 
     204          38 :     if (!m_oSRS.IsEmpty())
     205             :     {
     206           1 :         papszFileList = CSLAddString(papszFileList, osPrjFilename);
     207             :     }
     208             : 
     209          38 :     return papszFileList;
     210             : }
     211             : 
     212             : /************************************************************************/
     213             : /*                                Open()                                */
     214             : /************************************************************************/
     215             : 
     216          75 : GDALDataset *LCPDataset::Open(GDALOpenInfo *poOpenInfo)
     217             : 
     218             : {
     219             :     /* -------------------------------------------------------------------- */
     220             :     /*      Verify that this is a FARSITE LCP file                          */
     221             :     /* -------------------------------------------------------------------- */
     222          75 :     if (!Identify(poOpenInfo) || poOpenInfo->fpL == nullptr)
     223           0 :         return nullptr;
     224             : 
     225             :     /* -------------------------------------------------------------------- */
     226             :     /*      Confirm the requested access is supported.                      */
     227             :     /* -------------------------------------------------------------------- */
     228          75 :     if (poOpenInfo->eAccess == GA_Update)
     229             :     {
     230           0 :         ReportUpdateNotSupportedByDriver("LCP");
     231           0 :         return nullptr;
     232             :     }
     233             :     /* -------------------------------------------------------------------- */
     234             :     /*      Create a corresponding GDALDataset.                             */
     235             :     /* -------------------------------------------------------------------- */
     236         150 :     auto poDS = std::make_unique<LCPDataset>();
     237          75 :     std::swap(poDS->fpImage, poOpenInfo->fpL);
     238             : 
     239             :     /* -------------------------------------------------------------------- */
     240             :     /*      Read the header and extract some information.                   */
     241             :     /* -------------------------------------------------------------------- */
     242         150 :     if (VSIFSeekL(poDS->fpImage, 0, SEEK_SET) < 0 ||
     243          75 :         VSIFReadL(poDS->pachHeader, 1, LCP_HEADER_SIZE, poDS->fpImage) !=
     244             :             LCP_HEADER_SIZE)
     245             :     {
     246           0 :         CPLError(CE_Failure, CPLE_FileIO, "File too short");
     247           0 :         return nullptr;
     248             :     }
     249             : 
     250          75 :     const int nWidth = CPL_LSBSINT32PTR(poDS->pachHeader + 4164);
     251          75 :     const int nHeight = CPL_LSBSINT32PTR(poDS->pachHeader + 4168);
     252             : 
     253          75 :     poDS->nRasterXSize = nWidth;
     254          75 :     poDS->nRasterYSize = nHeight;
     255             : 
     256          75 :     if (!GDALCheckDatasetDimensions(poDS->nRasterXSize, poDS->nRasterYSize))
     257             :     {
     258           0 :         return nullptr;
     259             :     }
     260             : 
     261             :     // Crown fuels = canopy height, canopy base height, canopy bulk density.
     262             :     // 21 = have them, 20 = do not have them
     263             :     const bool bHaveCrownFuels =
     264          75 :         CPL_TO_BOOL(CPL_LSBSINT32PTR(poDS->pachHeader + 0) - 20);
     265             :     // Ground fuels = duff loading, coarse woody.
     266             :     const bool bHaveGroundFuels =
     267          75 :         CPL_TO_BOOL(CPL_LSBSINT32PTR(poDS->pachHeader + 4) - 20);
     268             : 
     269          75 :     int nBands = 0;
     270          75 :     if (bHaveCrownFuels)
     271             :     {
     272          69 :         if (bHaveGroundFuels)
     273          60 :             nBands = 10;
     274             :         else
     275           9 :             nBands = 8;
     276             :     }
     277             :     else
     278             :     {
     279           6 :         if (bHaveGroundFuels)
     280           3 :             nBands = 7;
     281             :         else
     282           3 :             nBands = 5;
     283             :     }
     284             : 
     285             :     // Add dataset-level metadata.
     286             : 
     287          75 :     int nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 8);
     288          75 :     char szTemp[32] = {'\0'};
     289          75 :     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     290          75 :     poDS->SetMetadataItem("LATITUDE", szTemp);
     291             : 
     292          75 :     nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 4204);
     293          75 :     if (nTemp == 0)
     294          74 :         poDS->SetMetadataItem("LINEAR_UNIT", "Meters");
     295          75 :     if (nTemp == 1)
     296           1 :         poDS->SetMetadataItem("LINEAR_UNIT", "Feet");
     297             : 
     298          75 :     poDS->pachHeader[LCP_HEADER_SIZE - 1] = '\0';
     299          75 :     poDS->SetMetadataItem("DESCRIPTION", poDS->pachHeader + 6804);
     300             : 
     301             :     /* -------------------------------------------------------------------- */
     302             :     /*      Create band information objects.                                */
     303             :     /* -------------------------------------------------------------------- */
     304             : 
     305          75 :     const int iPixelSize = nBands * 2;
     306             : 
     307          75 :     if (nWidth > INT_MAX / iPixelSize)
     308             :     {
     309           0 :         CPLError(CE_Failure, CPLE_AppDefined, "Int overflow occurred");
     310           0 :         return nullptr;
     311             :     }
     312             : 
     313         783 :     for (int iBand = 1; iBand <= nBands; iBand++)
     314             :     {
     315             :         auto poBand =
     316        1416 :             RawRasterBand::Create(poDS.get(), iBand, poDS->fpImage,
     317         708 :                                   LCP_HEADER_SIZE + ((iBand - 1) * 2),
     318             :                                   iPixelSize, iPixelSize * nWidth, GDT_Int16,
     319             :                                   RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN,
     320         708 :                                   RawRasterBand::OwnFP::NO);
     321         708 :         if (!poBand)
     322           0 :             return nullptr;
     323             : 
     324         708 :         switch (iBand)
     325             :         {
     326          75 :             case 1:
     327          75 :                 poBand->SetDescription("Elevation");
     328             : 
     329          75 :                 nTemp = CPL_LSBUINT16PTR(poDS->pachHeader + 4224);
     330          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     331          75 :                 poBand->SetMetadataItem("ELEVATION_UNIT", szTemp);
     332             : 
     333          75 :                 if (nTemp == 0)
     334          74 :                     poBand->SetMetadataItem("ELEVATION_UNIT_NAME", "Meters");
     335          75 :                 if (nTemp == 1)
     336           1 :                     poBand->SetMetadataItem("ELEVATION_UNIT_NAME", "Feet");
     337             : 
     338          75 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 44);
     339          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     340          75 :                 poBand->SetMetadataItem("ELEVATION_MIN", szTemp);
     341             : 
     342          75 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 48);
     343          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     344          75 :                 poBand->SetMetadataItem("ELEVATION_MAX", szTemp);
     345             : 
     346          75 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 52);
     347          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     348          75 :                 poBand->SetMetadataItem("ELEVATION_NUM_CLASSES", szTemp);
     349             : 
     350          75 :                 *(poDS->pachHeader + 4244 + 255) = '\0';
     351          75 :                 poBand->SetMetadataItem("ELEVATION_FILE",
     352          75 :                                         poDS->pachHeader + 4244);
     353             : 
     354          75 :                 break;
     355             : 
     356          75 :             case 2:
     357          75 :                 poBand->SetDescription("Slope");
     358             : 
     359          75 :                 nTemp = CPL_LSBUINT16PTR(poDS->pachHeader + 4226);
     360          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     361          75 :                 poBand->SetMetadataItem("SLOPE_UNIT", szTemp);
     362             : 
     363          75 :                 if (nTemp == 0)
     364          74 :                     poBand->SetMetadataItem("SLOPE_UNIT_NAME", "Degrees");
     365          75 :                 if (nTemp == 1)
     366           1 :                     poBand->SetMetadataItem("SLOPE_UNIT_NAME", "Percent");
     367             : 
     368          75 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 456);
     369          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     370          75 :                 poBand->SetMetadataItem("SLOPE_MIN", szTemp);
     371             : 
     372          75 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 460);
     373          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     374          75 :                 poBand->SetMetadataItem("SLOPE_MAX", szTemp);
     375             : 
     376          75 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 464);
     377          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     378          75 :                 poBand->SetMetadataItem("SLOPE_NUM_CLASSES", szTemp);
     379             : 
     380          75 :                 *(poDS->pachHeader + 4500 + 255) = '\0';
     381          75 :                 poBand->SetMetadataItem("SLOPE_FILE", poDS->pachHeader + 4500);
     382             : 
     383          75 :                 break;
     384             : 
     385          75 :             case 3:
     386          75 :                 poBand->SetDescription("Aspect");
     387             : 
     388          75 :                 nTemp = CPL_LSBUINT16PTR(poDS->pachHeader + 4228);
     389          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     390          75 :                 poBand->SetMetadataItem("ASPECT_UNIT", szTemp);
     391             : 
     392          75 :                 if (nTemp == 0)
     393           3 :                     poBand->SetMetadataItem("ASPECT_UNIT_NAME",
     394           3 :                                             "Grass categories");
     395          75 :                 if (nTemp == 1)
     396           1 :                     poBand->SetMetadataItem("ASPECT_UNIT_NAME",
     397           1 :                                             "Grass degrees");
     398          75 :                 if (nTemp == 2)
     399          71 :                     poBand->SetMetadataItem("ASPECT_UNIT_NAME",
     400          71 :                                             "Azimuth degrees");
     401             : 
     402          75 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 868);
     403          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     404          75 :                 poBand->SetMetadataItem("ASPECT_MIN", szTemp);
     405             : 
     406          75 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 872);
     407          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     408          75 :                 poBand->SetMetadataItem("ASPECT_MAX", szTemp);
     409             : 
     410          75 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 876);
     411          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     412          75 :                 poBand->SetMetadataItem("ASPECT_NUM_CLASSES", szTemp);
     413             : 
     414          75 :                 *(poDS->pachHeader + 4756 + 255) = '\0';
     415          75 :                 poBand->SetMetadataItem("ASPECT_FILE", poDS->pachHeader + 4756);
     416             : 
     417          75 :                 break;
     418             : 
     419          75 :             case 4:
     420             :             {
     421          75 :                 poBand->SetDescription("Fuel models");
     422             : 
     423          75 :                 nTemp = CPL_LSBUINT16PTR(poDS->pachHeader + 4230);
     424          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     425          75 :                 poBand->SetMetadataItem("FUEL_MODEL_OPTION", szTemp);
     426             : 
     427          75 :                 if (nTemp == 0)
     428          75 :                     poBand->SetMetadataItem(
     429             :                         "FUEL_MODEL_OPTION_DESC",
     430          75 :                         "no custom models AND no conversion file needed");
     431          75 :                 if (nTemp == 1)
     432           0 :                     poBand->SetMetadataItem(
     433             :                         "FUEL_MODEL_OPTION_DESC",
     434           0 :                         "custom models BUT no conversion file needed");
     435          75 :                 if (nTemp == 2)
     436           0 :                     poBand->SetMetadataItem(
     437             :                         "FUEL_MODEL_OPTION_DESC",
     438           0 :                         "no custom models BUT conversion file needed");
     439          75 :                 if (nTemp == 3)
     440           0 :                     poBand->SetMetadataItem(
     441             :                         "FUEL_MODEL_OPTION_DESC",
     442           0 :                         "custom models AND conversion file needed");
     443             : 
     444          75 :                 const int nMinFM = CPL_LSBSINT32PTR(poDS->pachHeader + 1280);
     445          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nMinFM);
     446          75 :                 poBand->SetMetadataItem("FUEL_MODEL_MIN", szTemp);
     447             : 
     448          75 :                 const int nMaxFM = CPL_LSBSINT32PTR(poDS->pachHeader + 1284);
     449          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nMaxFM);
     450          75 :                 poBand->SetMetadataItem("FUEL_MODEL_MAX", szTemp);
     451             : 
     452          75 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 1288);
     453          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     454          75 :                 poBand->SetMetadataItem("FUEL_MODEL_NUM_CLASSES", szTemp);
     455             : 
     456         150 :                 std::string osValues;
     457          75 :                 if (nTemp > 0 && nTemp <= 100)
     458             :                 {
     459         267 :                     for (int i = 0; i <= nTemp; i++)
     460             :                     {
     461         192 :                         const int nTemp2 = CPL_LSBSINT32PTR(poDS->pachHeader +
     462             :                                                             (1292 + (i * 4)));
     463         192 :                         if (nTemp2 >= nMinFM && nTemp2 <= nMaxFM)
     464             :                         {
     465         100 :                             snprintf(szTemp, sizeof(szTemp), "%d", nTemp2);
     466         100 :                             if (!osValues.empty())
     467          27 :                                 osValues += ',';
     468         100 :                             osValues += szTemp;
     469             :                         }
     470             :                     }
     471             :                 }
     472          75 :                 poBand->SetMetadataItem("FUEL_MODEL_VALUES", osValues.c_str());
     473             : 
     474          75 :                 *(poDS->pachHeader + 5012 + 255) = '\0';
     475          75 :                 poBand->SetMetadataItem("FUEL_MODEL_FILE",
     476          75 :                                         poDS->pachHeader + 5012);
     477             : 
     478          75 :                 break;
     479             :             }
     480          75 :             case 5:
     481          75 :                 poBand->SetDescription("Canopy cover");
     482             : 
     483          75 :                 nTemp = CPL_LSBUINT16PTR(poDS->pachHeader + 4232);
     484          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     485          75 :                 poBand->SetMetadataItem("CANOPY_COV_UNIT", szTemp);
     486             : 
     487          75 :                 if (nTemp == 0)
     488           4 :                     poBand->SetMetadataItem("CANOPY_COV_UNIT_NAME",
     489           4 :                                             "Categories (0-4)");
     490          75 :                 if (nTemp == 1)
     491          71 :                     poBand->SetMetadataItem("CANOPY_COV_UNIT_NAME", "Percent");
     492             : 
     493          75 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 1692);
     494          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     495          75 :                 poBand->SetMetadataItem("CANOPY_COV_MIN", szTemp);
     496             : 
     497          75 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 1696);
     498          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     499          75 :                 poBand->SetMetadataItem("CANOPY_COV_MAX", szTemp);
     500             : 
     501          75 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 1700);
     502          75 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     503          75 :                 poBand->SetMetadataItem("CANOPY_COV_NUM_CLASSES", szTemp);
     504             : 
     505          75 :                 *(poDS->pachHeader + 5268 + 255) = '\0';
     506          75 :                 poBand->SetMetadataItem("CANOPY_COV_FILE",
     507          75 :                                         poDS->pachHeader + 5268);
     508             : 
     509          75 :                 break;
     510             : 
     511          72 :             case 6:
     512          72 :                 if (bHaveCrownFuels)
     513             :                 {
     514          69 :                     poBand->SetDescription("Canopy height");
     515             : 
     516          69 :                     nTemp = CPL_LSBUINT16PTR(poDS->pachHeader + 4234);
     517          69 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     518          69 :                     poBand->SetMetadataItem("CANOPY_HT_UNIT", szTemp);
     519             : 
     520          69 :                     if (nTemp == 1)
     521           3 :                         poBand->SetMetadataItem("CANOPY_HT_UNIT_NAME",
     522           3 :                                                 "Meters");
     523          69 :                     if (nTemp == 2)
     524           3 :                         poBand->SetMetadataItem("CANOPY_HT_UNIT_NAME", "Feet");
     525          69 :                     if (nTemp == 3)
     526          62 :                         poBand->SetMetadataItem("CANOPY_HT_UNIT_NAME",
     527          62 :                                                 "Meters x 10");
     528          69 :                     if (nTemp == 4)
     529           1 :                         poBand->SetMetadataItem("CANOPY_HT_UNIT_NAME",
     530           1 :                                                 "Feet x 10");
     531             : 
     532          69 :                     nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 2104);
     533          69 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     534          69 :                     poBand->SetMetadataItem("CANOPY_HT_MIN", szTemp);
     535             : 
     536          69 :                     nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 2108);
     537          69 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     538          69 :                     poBand->SetMetadataItem("CANOPY_HT_MAX", szTemp);
     539             : 
     540          69 :                     nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 2112);
     541          69 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     542          69 :                     poBand->SetMetadataItem("CANOPY_HT_NUM_CLASSES", szTemp);
     543             : 
     544          69 :                     *(poDS->pachHeader + 5524 + 255) = '\0';
     545          69 :                     poBand->SetMetadataItem("CANOPY_HT_FILE",
     546          69 :                                             poDS->pachHeader + 5524);
     547             :                 }
     548             :                 else
     549             :                 {
     550           3 :                     poBand->SetDescription("Duff");
     551             : 
     552           3 :                     nTemp = CPL_LSBUINT16PTR(poDS->pachHeader + 4240);
     553           3 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     554           3 :                     poBand->SetMetadataItem("DUFF_UNIT", szTemp);
     555             : 
     556           3 :                     if (nTemp == 1)
     557           3 :                         poBand->SetMetadataItem("DUFF_UNIT_NAME", "Mg/ha");
     558           3 :                     if (nTemp == 2)
     559           0 :                         poBand->SetMetadataItem("DUFF_UNIT_NAME", "t/ac");
     560             : 
     561           3 :                     nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 3340);
     562           3 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     563           3 :                     poBand->SetMetadataItem("DUFF_MIN", szTemp);
     564             : 
     565           3 :                     nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 3344);
     566           3 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     567           3 :                     poBand->SetMetadataItem("DUFF_MAX", szTemp);
     568             : 
     569           3 :                     nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 3348);
     570           3 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     571           3 :                     poBand->SetMetadataItem("DUFF_NUM_CLASSES", szTemp);
     572             : 
     573           3 :                     *(poDS->pachHeader + 6292 + 255) = '\0';
     574           3 :                     poBand->SetMetadataItem("DUFF_FILE",
     575           3 :                                             poDS->pachHeader + 6292);
     576             :                 }
     577          72 :                 break;
     578             : 
     579          72 :             case 7:
     580          72 :                 if (bHaveCrownFuels)
     581             :                 {
     582          69 :                     poBand->SetDescription("Canopy base height");
     583             : 
     584          69 :                     nTemp = CPL_LSBUINT16PTR(poDS->pachHeader + 4236);
     585          69 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     586          69 :                     poBand->SetMetadataItem("CBH_UNIT", szTemp);
     587             : 
     588          69 :                     if (nTemp == 1)
     589           3 :                         poBand->SetMetadataItem("CBH_UNIT_NAME", "Meters");
     590          69 :                     if (nTemp == 2)
     591           3 :                         poBand->SetMetadataItem("CBH_UNIT_NAME", "Feet");
     592          69 :                     if (nTemp == 3)
     593          62 :                         poBand->SetMetadataItem("CBH_UNIT_NAME", "Meters x 10");
     594          69 :                     if (nTemp == 4)
     595           1 :                         poBand->SetMetadataItem("CBH_UNIT_NAME", "Feet x 10");
     596             : 
     597          69 :                     nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 2516);
     598          69 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     599          69 :                     poBand->SetMetadataItem("CBH_MIN", szTemp);
     600             : 
     601          69 :                     nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 2520);
     602          69 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     603          69 :                     poBand->SetMetadataItem("CBH_MAX", szTemp);
     604             : 
     605          69 :                     nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 2524);
     606          69 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     607          69 :                     poBand->SetMetadataItem("CBH_NUM_CLASSES", szTemp);
     608             : 
     609          69 :                     *(poDS->pachHeader + 5780 + 255) = '\0';
     610          69 :                     poBand->SetMetadataItem("CBH_FILE",
     611          69 :                                             poDS->pachHeader + 5780);
     612             :                 }
     613             :                 else
     614             :                 {
     615           3 :                     poBand->SetDescription("Coarse woody debris");
     616             : 
     617           3 :                     nTemp = CPL_LSBUINT16PTR(poDS->pachHeader + 4242);
     618           3 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     619           3 :                     poBand->SetMetadataItem("CWD_OPTION", szTemp);
     620             : 
     621             :                     // if ( nTemp == 1 )
     622             :                     //    poBand->SetMetadataItem( "CWD_UNIT_DESC", "?" );
     623             :                     // if ( nTemp == 2 )
     624             :                     //    poBand->SetMetadataItem( "CWD_UNIT_DESC", "?" );
     625             : 
     626           3 :                     nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 3752);
     627           3 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     628           3 :                     poBand->SetMetadataItem("CWD_MIN", szTemp);
     629             : 
     630           3 :                     nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 3756);
     631           3 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     632           3 :                     poBand->SetMetadataItem("CWD_MAX", szTemp);
     633             : 
     634           3 :                     nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 3760);
     635           3 :                     snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     636           3 :                     poBand->SetMetadataItem("CWD_NUM_CLASSES", szTemp);
     637             : 
     638           3 :                     *(poDS->pachHeader + 6548 + 255) = '\0';
     639           3 :                     poBand->SetMetadataItem("CWD_FILE",
     640           3 :                                             poDS->pachHeader + 6548);
     641             :                 }
     642          72 :                 break;
     643             : 
     644          69 :             case 8:
     645          69 :                 poBand->SetDescription("Canopy bulk density");
     646             : 
     647          69 :                 nTemp = CPL_LSBUINT16PTR(poDS->pachHeader + 4238);
     648          69 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     649          69 :                 poBand->SetMetadataItem("CBD_UNIT", szTemp);
     650             : 
     651          69 :                 if (nTemp == 1)
     652           3 :                     poBand->SetMetadataItem("CBD_UNIT_NAME", "kg/m^3");
     653          69 :                 if (nTemp == 2)
     654           3 :                     poBand->SetMetadataItem("CBD_UNIT_NAME", "lb/ft^3");
     655          69 :                 if (nTemp == 3)
     656          62 :                     poBand->SetMetadataItem("CBD_UNIT_NAME", "kg/m^3 x 100");
     657          69 :                 if (nTemp == 4)
     658           1 :                     poBand->SetMetadataItem("CBD_UNIT_NAME", "lb/ft^3 x 1000");
     659             : 
     660          69 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 2928);
     661          69 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     662          69 :                 poBand->SetMetadataItem("CBD_MIN", szTemp);
     663             : 
     664          69 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 2932);
     665          69 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     666          69 :                 poBand->SetMetadataItem("CBD_MAX", szTemp);
     667             : 
     668          69 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 2936);
     669          69 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     670          69 :                 poBand->SetMetadataItem("CBD_NUM_CLASSES", szTemp);
     671             : 
     672          69 :                 *(poDS->pachHeader + 6036 + 255) = '\0';
     673          69 :                 poBand->SetMetadataItem("CBD_FILE", poDS->pachHeader + 6036);
     674             : 
     675          69 :                 break;
     676             : 
     677          60 :             case 9:
     678          60 :                 poBand->SetDescription("Duff");
     679             : 
     680          60 :                 nTemp = CPL_LSBUINT16PTR(poDS->pachHeader + 4240);
     681          60 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     682          60 :                 poBand->SetMetadataItem("DUFF_UNIT", szTemp);
     683             : 
     684          60 :                 if (nTemp == 1)
     685          59 :                     poBand->SetMetadataItem("DUFF_UNIT_NAME", "Mg/ha");
     686          60 :                 if (nTemp == 2)
     687           1 :                     poBand->SetMetadataItem("DUFF_UNIT_NAME", "t/ac");
     688             : 
     689          60 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 3340);
     690          60 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     691          60 :                 poBand->SetMetadataItem("DUFF_MIN", szTemp);
     692             : 
     693          60 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 3344);
     694          60 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     695          60 :                 poBand->SetMetadataItem("DUFF_MAX", szTemp);
     696             : 
     697          60 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 3348);
     698          60 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     699          60 :                 poBand->SetMetadataItem("DUFF_NUM_CLASSES", szTemp);
     700             : 
     701          60 :                 *(poDS->pachHeader + 6292 + 255) = '\0';
     702          60 :                 poBand->SetMetadataItem("DUFF_FILE", poDS->pachHeader + 6292);
     703             : 
     704          60 :                 break;
     705             : 
     706          60 :             case 10:
     707          60 :                 poBand->SetDescription("Coarse woody debris");
     708             : 
     709          60 :                 nTemp = CPL_LSBUINT16PTR(poDS->pachHeader + 4242);
     710          60 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     711          60 :                 poBand->SetMetadataItem("CWD_OPTION", szTemp);
     712             : 
     713             :                 // if ( nTemp == 1 )
     714             :                 //    poBand->SetMetadataItem( "CWD_UNIT_DESC", "?" );
     715             :                 // if ( nTemp == 2 )
     716             :                 //    poBand->SetMetadataItem( "CWD_UNIT_DESC", "?" );
     717             : 
     718          60 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 3752);
     719          60 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     720          60 :                 poBand->SetMetadataItem("CWD_MIN", szTemp);
     721             : 
     722          60 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 3756);
     723          60 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     724          60 :                 poBand->SetMetadataItem("CWD_MAX", szTemp);
     725             : 
     726          60 :                 nTemp = CPL_LSBSINT32PTR(poDS->pachHeader + 3760);
     727          60 :                 snprintf(szTemp, sizeof(szTemp), "%d", nTemp);
     728          60 :                 poBand->SetMetadataItem("CWD_NUM_CLASSES", szTemp);
     729             : 
     730          60 :                 *(poDS->pachHeader + 6548 + 255) = '\0';
     731          60 :                 poBand->SetMetadataItem("CWD_FILE", poDS->pachHeader + 6548);
     732             : 
     733          60 :                 break;
     734             :         }
     735             : 
     736         708 :         poDS->SetBand(iBand, std::move(poBand));
     737             :     }
     738             : 
     739             :     /* -------------------------------------------------------------------- */
     740             :     /*      Try to read projection file.                                    */
     741             :     /* -------------------------------------------------------------------- */
     742             :     char *const pszDirname =
     743          75 :         CPLStrdup(CPLGetPathSafe(poOpenInfo->pszFilename).c_str());
     744             :     char *const pszBasename =
     745          75 :         CPLStrdup(CPLGetBasenameSafe(poOpenInfo->pszFilename).c_str());
     746             : 
     747          75 :     poDS->osPrjFilename = CPLFormFilenameSafe(pszDirname, pszBasename, "prj");
     748             :     VSIStatBufL sStatBuf;
     749          75 :     int nRet = VSIStatL(poDS->osPrjFilename, &sStatBuf);
     750             : 
     751          75 :     if (nRet != 0 && VSIIsCaseSensitiveFS(poDS->osPrjFilename))
     752             :     {
     753          72 :         poDS->osPrjFilename =
     754         144 :             CPLFormFilenameSafe(pszDirname, pszBasename, "PRJ");
     755          72 :         nRet = VSIStatL(poDS->osPrjFilename, &sStatBuf);
     756             :     }
     757             : 
     758          75 :     if (nRet == 0)
     759             :     {
     760           3 :         char **papszPrj = CSLLoad(poDS->osPrjFilename);
     761             : 
     762           3 :         CPLDebug("LCP", "Loaded SRS from %s", poDS->osPrjFilename.c_str());
     763             : 
     764           3 :         poDS->m_oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
     765           3 :         if (poDS->m_oSRS.importFromESRI(papszPrj) != OGRERR_NONE)
     766             :         {
     767           0 :             poDS->m_oSRS.Clear();
     768             :         }
     769             : 
     770           3 :         CSLDestroy(papszPrj);
     771             :     }
     772             : 
     773          75 :     CPLFree(pszDirname);
     774          75 :     CPLFree(pszBasename);
     775             : 
     776             :     /* -------------------------------------------------------------------- */
     777             :     /*      Initialize any PAM information.                                 */
     778             :     /* -------------------------------------------------------------------- */
     779          75 :     poDS->SetDescription(poOpenInfo->pszFilename);
     780          75 :     poDS->TryLoadXML();
     781             : 
     782             :     /* -------------------------------------------------------------------- */
     783             :     /*      Check for external overviews.                                   */
     784             :     /* -------------------------------------------------------------------- */
     785         150 :     poDS->oOvManager.Initialize(poDS.get(), poOpenInfo->pszFilename,
     786          75 :                                 poOpenInfo->GetSiblingFiles());
     787             : 
     788          75 :     return poDS.release();
     789             : }
     790             : 
     791             : /************************************************************************/
     792             : /*                          ClassifyBandData()                          */
     793             : /*  Classify a band and store 99 or fewer unique values.  If there are  */
     794             : /*  more than 99 unique values, then set pnNumClasses to -1 as a flag   */
     795             : /*  that represents this.  These are legacy values in the header, and   */
     796             : /*  while we should never deprecate them, we could possibly not         */
     797             : /*  calculate them by default.                                          */
     798             : /************************************************************************/
     799             : 
     800         320 : CPLErr LCPDataset::ClassifyBandData(GDALRasterBand *poBand, GInt32 &nNumClasses,
     801             :                                     GInt32 *panClasses)
     802             : {
     803         320 :     CPLAssert(poBand);
     804         320 :     CPLAssert(panClasses);
     805             : 
     806         320 :     const int nXSize = poBand->GetXSize();
     807         320 :     const int nYSize = poBand->GetYSize();
     808             : 
     809             :     GInt16 *panValues =
     810         320 :         static_cast<GInt16 *>(CPLMalloc(sizeof(GInt16) * nXSize));
     811         320 :     constexpr int MIN_VAL = std::numeric_limits<GInt16>::min();
     812         320 :     constexpr int MAX_VAL = std::numeric_limits<GInt16>::max();
     813         320 :     constexpr int RANGE_VAL = MAX_VAL - MIN_VAL + 1;
     814         320 :     GByte *pabyFlags = static_cast<GByte *>(CPLCalloc(1, RANGE_VAL));
     815             : 
     816             :     /* so that values in [-32768,32767] are mapped to [0,65535] in pabyFlags */
     817         320 :     constexpr int OFFSET = -MIN_VAL;
     818             : 
     819         320 :     int nFound = 0;
     820         320 :     bool bTooMany = false;
     821         320 :     CPLErr eErr = CE_None;
     822        3770 :     for (int iLine = 0; iLine < nYSize; iLine++)
     823             :     {
     824        3451 :         eErr = poBand->RasterIO(GF_Read, 0, iLine, nXSize, 1, panValues, nXSize,
     825             :                                 1, GDT_Int16, 0, 0, nullptr);
     826        3451 :         if (eErr != CE_None)
     827           0 :             break;
     828       37531 :         for (int iPixel = 0; iPixel < nXSize; iPixel++)
     829             :         {
     830       34081 :             if (panValues[iPixel] == -9999)
     831             :             {
     832           0 :                 continue;
     833             :             }
     834       34081 :             if (nFound == LCP_MAX_CLASSES)
     835             :             {
     836           1 :                 CPLDebug("LCP",
     837             :                          "Found more that %d unique values in "
     838             :                          "band %d.  Not 'classifying' the data.",
     839             :                          LCP_MAX_CLASSES - 1, poBand->GetBand());
     840           1 :                 nFound = -1;
     841           1 :                 bTooMany = true;
     842           1 :                 break;
     843             :             }
     844       34080 :             if (pabyFlags[panValues[iPixel] + OFFSET] == 0)
     845             :             {
     846         575 :                 pabyFlags[panValues[iPixel] + OFFSET] = 1;
     847         575 :                 nFound++;
     848             :             }
     849             :         }
     850        3451 :         if (bTooMany)
     851           1 :             break;
     852             :     }
     853         320 :     if (!bTooMany)
     854             :     {
     855             :         // The classes are always padded with a leading 0.  This was for aligning
     856             :         // offsets, or making it a 1-based array instead of 0-based.
     857         319 :         panClasses[0] = 0;
     858    20906300 :         for (int j = 0, nIndex = 1; j < RANGE_VAL; j++)
     859             :         {
     860    20906000 :             if (pabyFlags[j] == 1)
     861             :             {
     862         475 :                 panClasses[nIndex++] = j;
     863             :             }
     864             :         }
     865             :     }
     866         320 :     nNumClasses = nFound;
     867         320 :     CPLFree(pabyFlags);
     868         320 :     CPLFree(panValues);
     869             : 
     870         320 :     return eErr;
     871             : }
     872             : 
     873             : /************************************************************************/
     874             : /*                          CreateCopy()                                */
     875             : /************************************************************************/
     876             : 
     877          68 : GDALDataset *LCPDataset::CreateCopy(const char *pszFilename,
     878             :                                     GDALDataset *poSrcDS, int bStrict,
     879             :                                     char **papszOptions,
     880             :                                     GDALProgressFunc pfnProgress,
     881             :                                     void *pProgressData)
     882             : 
     883             : {
     884             :     /* -------------------------------------------------------------------- */
     885             :     /*      Verify input options.                                           */
     886             :     /* -------------------------------------------------------------------- */
     887          68 :     const int nBands = poSrcDS->GetRasterCount();
     888          68 :     if (nBands != 5 && nBands != 7 && nBands != 8 && nBands != 10)
     889             :     {
     890          25 :         CPLError(CE_Failure, CPLE_NotSupported,
     891             :                  "LCP driver doesn't support %d bands.  Must be 5, 7, 8 "
     892             :                  "or 10 bands.",
     893             :                  nBands);
     894          25 :         return nullptr;
     895             :     }
     896             : 
     897          43 :     GDALDataType eType = poSrcDS->GetRasterBand(1)->GetRasterDataType();
     898          43 :     if (eType != GDT_Int16 && bStrict)
     899             :     {
     900           1 :         CPLError(CE_Failure, CPLE_AppDefined,
     901             :                  "LCP only supports 16-bit signed integer data types.");
     902           1 :         return nullptr;
     903             :     }
     904          42 :     else if (eType != GDT_Int16)
     905             :     {
     906           0 :         CPLError(CE_Warning, CPLE_AppDefined,
     907             :                  "Setting data type to 16-bit integer.");
     908             :     }
     909             : 
     910             :     /* -------------------------------------------------------------------- */
     911             :     /*      What schema do we have (ground/crown fuels)                     */
     912             :     /* -------------------------------------------------------------------- */
     913             : 
     914          42 :     const bool bHaveCrownFuels = nBands == 8 || nBands == 10;
     915          42 :     const bool bHaveGroundFuels = nBands == 7 || nBands == 10;
     916             : 
     917             :     // Since units are 'configurable', we should check for user
     918             :     // defined units.  This is a bit cumbersome, but the user should
     919             :     // be allowed to specify none to get default units/options.  Use
     920             :     // default units every chance we get.
     921          42 :     GInt16 panMetadata[LCP_MAX_BANDS] = {
     922             :         0,  // 0 ELEVATION_UNIT
     923             :         0,  // 1 SLOPE_UNIT
     924             :         2,  // 2 ASPECT_UNIT
     925             :         0,  // 3 FUEL_MODEL_OPTION
     926             :         1,  // 4 CANOPY_COV_UNIT
     927             :         3,  // 5 CANOPY_HT_UNIT
     928             :         3,  // 6 CBH_UNIT
     929             :         3,  // 7 CBD_UNIT
     930             :         1,  // 8 DUFF_UNIT
     931             :         0,  // 9 CWD_OPTION
     932             :     };
     933             : 
     934             :     // Check the units/options for user overrides.
     935             :     const char *pszTemp =
     936          42 :         CSLFetchNameValueDef(papszOptions, "ELEVATION_UNIT", "METERS");
     937          42 :     if (STARTS_WITH_CI(pszTemp, "METER"))
     938             :     {
     939          40 :         panMetadata[0] = 0;
     940             :     }
     941           2 :     else if (EQUAL(pszTemp, "FEET") || EQUAL(pszTemp, "FOOT"))
     942             :     {
     943           1 :         panMetadata[0] = 1;
     944             :     }
     945             :     else
     946             :     {
     947           1 :         CPLError(CE_Failure, CPLE_AppDefined,
     948             :                  "Invalid value (%s) for ELEVATION_UNIT.", pszTemp);
     949           1 :         return nullptr;
     950             :     }
     951             : 
     952          41 :     pszTemp = CSLFetchNameValueDef(papszOptions, "SLOPE_UNIT", "DEGREES");
     953          41 :     if (EQUAL(pszTemp, "DEGREES"))
     954             :     {
     955          39 :         panMetadata[1] = 0;
     956             :     }
     957           2 :     else if (EQUAL(pszTemp, "PERCENT"))
     958             :     {
     959           1 :         panMetadata[1] = 1;
     960             :     }
     961             :     else
     962             :     {
     963           1 :         CPLError(CE_Failure, CPLE_AppDefined,
     964             :                  "Invalid value (%s) for SLOPE_UNIT.", pszTemp);
     965           1 :         return nullptr;
     966             :     }
     967             : 
     968             :     pszTemp =
     969          40 :         CSLFetchNameValueDef(papszOptions, "ASPECT_UNIT", "AZIMUTH_DEGREES");
     970          40 :     if (EQUAL(pszTemp, "GRASS_CATEGORIES"))
     971             :     {
     972           1 :         panMetadata[2] = 0;
     973             :     }
     974          39 :     else if (EQUAL(pszTemp, "GRASS_DEGREES"))
     975             :     {
     976           1 :         panMetadata[2] = 1;
     977             :     }
     978          38 :     else if (EQUAL(pszTemp, "AZIMUTH_DEGREES"))
     979             :     {
     980          37 :         panMetadata[2] = 2;
     981             :     }
     982             :     else
     983             :     {
     984           1 :         CPLError(CE_Failure, CPLE_AppDefined,
     985             :                  "Invalid value (%s) for ASPECT_UNIT.", pszTemp);
     986           1 :         return nullptr;
     987             :     }
     988             : 
     989          39 :     pszTemp = CSLFetchNameValueDef(papszOptions, "FUEL_MODEL_OPTION",
     990             :                                    "NO_CUSTOM_AND_NO_FILE");
     991          39 :     if (EQUAL(pszTemp, "NO_CUSTOM_AND_NO_FILE"))
     992             :     {
     993          38 :         panMetadata[3] = 0;
     994             :     }
     995           1 :     else if (EQUAL(pszTemp, "CUSTOM_AND_NO_FILE"))
     996             :     {
     997           0 :         panMetadata[3] = 1;
     998             :     }
     999           1 :     else if (EQUAL(pszTemp, "NO_CUSTOM_AND_FILE"))
    1000             :     {
    1001           0 :         panMetadata[3] = 2;
    1002             :     }
    1003           1 :     else if (EQUAL(pszTemp, "CUSTOM_AND_FILE"))
    1004             :     {
    1005           0 :         panMetadata[3] = 3;
    1006             :     }
    1007             :     else
    1008             :     {
    1009           1 :         CPLError(CE_Failure, CPLE_AppDefined,
    1010             :                  "Invalid value (%s) for FUEL_MODEL_OPTION.", pszTemp);
    1011           1 :         return nullptr;
    1012             :     }
    1013             : 
    1014          38 :     pszTemp = CSLFetchNameValueDef(papszOptions, "CANOPY_COV_UNIT", "PERCENT");
    1015          38 :     if (EQUAL(pszTemp, "CATEGORIES"))
    1016             :     {
    1017           1 :         panMetadata[4] = 0;
    1018             :     }
    1019          37 :     else if (EQUAL(pszTemp, "PERCENT"))
    1020             :     {
    1021          36 :         panMetadata[4] = 1;
    1022             :     }
    1023             :     else
    1024             :     {
    1025           1 :         CPLError(CE_Failure, CPLE_AppDefined,
    1026             :                  "Invalid value (%s) for CANOPY_COV_UNIT.", pszTemp);
    1027           1 :         return nullptr;
    1028             :     }
    1029             : 
    1030          37 :     if (bHaveCrownFuels)
    1031             :     {
    1032             :         pszTemp =
    1033          35 :             CSLFetchNameValueDef(papszOptions, "CANOPY_HT_UNIT", "METERS_X_10");
    1034          35 :         if (EQUAL(pszTemp, "METERS") || EQUAL(pszTemp, "METER"))
    1035             :         {
    1036           1 :             panMetadata[5] = 1;
    1037             :         }
    1038          34 :         else if (EQUAL(pszTemp, "FEET") || EQUAL(pszTemp, "FOOT"))
    1039             :         {
    1040           1 :             panMetadata[5] = 2;
    1041             :         }
    1042          33 :         else if (EQUAL(pszTemp, "METERS_X_10") || EQUAL(pszTemp, "METER_X_10"))
    1043             :         {
    1044          31 :             panMetadata[5] = 3;
    1045             :         }
    1046           2 :         else if (EQUAL(pszTemp, "FEET_X_10") || EQUAL(pszTemp, "FOOT_X_10"))
    1047             :         {
    1048           1 :             panMetadata[5] = 4;
    1049             :         }
    1050             :         else
    1051             :         {
    1052           1 :             CPLError(CE_Failure, CPLE_AppDefined,
    1053             :                      "Invalid value (%s) for CANOPY_HT_UNIT.", pszTemp);
    1054           1 :             return nullptr;
    1055             :         }
    1056             : 
    1057          34 :         pszTemp = CSLFetchNameValueDef(papszOptions, "CBH_UNIT", "METERS_X_10");
    1058          34 :         if (EQUAL(pszTemp, "METERS") || EQUAL(pszTemp, "METER"))
    1059             :         {
    1060           1 :             panMetadata[6] = 1;
    1061             :         }
    1062          33 :         else if (EQUAL(pszTemp, "FEET") || EQUAL(pszTemp, "FOOT"))
    1063             :         {
    1064           1 :             panMetadata[6] = 2;
    1065             :         }
    1066          32 :         else if (EQUAL(pszTemp, "METERS_X_10") || EQUAL(pszTemp, "METER_X_10"))
    1067             :         {
    1068          30 :             panMetadata[6] = 3;
    1069             :         }
    1070           2 :         else if (EQUAL(pszTemp, "FEET_X_10") || EQUAL(pszTemp, "FOOT_X_10"))
    1071             :         {
    1072           1 :             panMetadata[6] = 4;
    1073             :         }
    1074             :         else
    1075             :         {
    1076           1 :             CPLError(CE_Failure, CPLE_AppDefined,
    1077             :                      "Invalid value (%s) for CBH_UNIT.", pszTemp);
    1078           1 :             return nullptr;
    1079             :         }
    1080             : 
    1081          33 :         pszTemp = CSLFetchNameValueDef(papszOptions, "CBD_UNIT",
    1082             :                                        "KG_PER_CUBIC_METER_X_100");
    1083          33 :         if (EQUAL(pszTemp, "KG_PER_CUBIC_METER"))
    1084             :         {
    1085           1 :             panMetadata[7] = 1;
    1086             :         }
    1087          32 :         else if (EQUAL(pszTemp, "POUND_PER_CUBIC_FOOT"))
    1088             :         {
    1089           1 :             panMetadata[7] = 2;
    1090             :         }
    1091          31 :         else if (EQUAL(pszTemp, "KG_PER_CUBIC_METER_X_100"))
    1092             :         {
    1093          29 :             panMetadata[7] = 3;
    1094             :         }
    1095           2 :         else if (EQUAL(pszTemp, "POUND_PER_CUBIC_FOOT_X_1000"))
    1096             :         {
    1097           1 :             panMetadata[7] = 4;
    1098             :         }
    1099             :         else
    1100             :         {
    1101           1 :             CPLError(CE_Failure, CPLE_AppDefined,
    1102             :                      "Invalid value (%s) for CBD_UNIT.", pszTemp);
    1103           1 :             return nullptr;
    1104             :         }
    1105             :     }
    1106             : 
    1107          34 :     if (bHaveGroundFuels)
    1108             :     {
    1109          32 :         pszTemp = CSLFetchNameValueDef(papszOptions, "DUFF_UNIT",
    1110             :                                        "MG_PER_HECTARE_X_10");
    1111          32 :         if (EQUAL(pszTemp, "MG_PER_HECTARE_X_10"))
    1112             :         {
    1113          30 :             panMetadata[8] = 1;
    1114             :         }
    1115           2 :         else if (EQUAL(pszTemp, "TONS_PER_ACRE_X_10"))
    1116             :         {
    1117           1 :             panMetadata[8] = 2;
    1118             :         }
    1119             :         else
    1120             :         {
    1121           1 :             CPLError(CE_Failure, CPLE_AppDefined,
    1122             :                      "Invalid value (%s) for DUFF_UNIT.", pszTemp);
    1123           1 :             return nullptr;
    1124             :         }
    1125             : 
    1126          31 :         panMetadata[9] = 1;
    1127             :     }
    1128             : 
    1129             :     // Calculate the stats for each band.  The binary file carries along
    1130             :     // these metadata for display purposes(?).
    1131          33 :     bool bCalculateStats = CPLFetchBool(papszOptions, "CALCULATE_STATS", true);
    1132             : 
    1133             :     const bool bClassifyData =
    1134          33 :         CPLFetchBool(papszOptions, "CLASSIFY_DATA", true);
    1135             : 
    1136             :     // We should have stats if we classify, we'll get them anyway.
    1137          33 :     if (bClassifyData && !bCalculateStats)
    1138             :     {
    1139           0 :         CPLError(CE_Warning, CPLE_AppDefined,
    1140             :                  "Ignoring request to not calculate statistics, "
    1141             :                  "because CLASSIFY_DATA was set to ON");
    1142           0 :         bCalculateStats = true;
    1143             :     }
    1144             : 
    1145          33 :     pszTemp = CSLFetchNameValueDef(papszOptions, "LINEAR_UNIT", "SET_FROM_SRS");
    1146          33 :     int nLinearUnits = 0;
    1147          33 :     bool bSetLinearUnits = false;
    1148          33 :     if (EQUAL(pszTemp, "SET_FROM_SRS"))
    1149             :     {
    1150           0 :         bSetLinearUnits = true;
    1151             :     }
    1152          33 :     else if (STARTS_WITH_CI(pszTemp, "METER"))
    1153             :     {
    1154          32 :         nLinearUnits = 0;
    1155             :     }
    1156           1 :     else if (EQUAL(pszTemp, "FOOT") || EQUAL(pszTemp, "FEET"))
    1157             :     {
    1158           1 :         nLinearUnits = 1;
    1159             :     }
    1160           0 :     else if (STARTS_WITH_CI(pszTemp, "KILOMETER"))
    1161             :     {
    1162           0 :         nLinearUnits = 2;
    1163             :     }
    1164          33 :     bool bCalculateLatitude = true;
    1165          33 :     int nLatitude = 0;
    1166          33 :     if (CSLFetchNameValue(papszOptions, "LATITUDE") != nullptr)
    1167             :     {
    1168          33 :         bCalculateLatitude = false;
    1169          33 :         nLatitude = atoi(CSLFetchNameValue(papszOptions, "LATITUDE"));
    1170          33 :         if (nLatitude > 90 || nLatitude < -90)
    1171             :         {
    1172           0 :             CPLError(CE_Failure, CPLE_OpenFailed,
    1173             :                      "Invalid value (%d) for LATITUDE.", nLatitude);
    1174           0 :             return nullptr;
    1175             :         }
    1176             :     }
    1177             : 
    1178             :     // If no latitude is supplied, attempt to extract the central latitude
    1179             :     // from the image.  It must be set either manually or here, otherwise
    1180             :     // we fail.
    1181          33 :     GDALGeoTransform srcGT;
    1182          33 :     poSrcDS->GetGeoTransform(srcGT);
    1183          33 :     const OGRSpatialReference *poSrcSRS = poSrcDS->GetSpatialRef();
    1184          33 :     double dfLongitude = 0.0;
    1185          33 :     double dfLatitude = 0.0;
    1186             : 
    1187          33 :     const int nYSize = poSrcDS->GetRasterYSize();
    1188             : 
    1189          33 :     if (!bCalculateLatitude)
    1190             :     {
    1191          33 :         dfLatitude = nLatitude;
    1192             :     }
    1193           0 :     else if (poSrcSRS)
    1194             :     {
    1195           0 :         OGRSpatialReference oDstSRS;
    1196           0 :         oDstSRS.importFromEPSG(4269);
    1197           0 :         oDstSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
    1198             :         OGRCoordinateTransformation *poCT =
    1199           0 :             reinterpret_cast<OGRCoordinateTransformation *>(
    1200             :                 OGRCreateCoordinateTransformation(poSrcSRS, &oDstSRS));
    1201           0 :         if (poCT != nullptr)
    1202             :         {
    1203           0 :             dfLatitude = srcGT[3] + srcGT[5] * nYSize / 2;
    1204             :             const int nErr =
    1205           0 :                 static_cast<int>(poCT->Transform(1, &dfLongitude, &dfLatitude));
    1206           0 :             if (!nErr)
    1207             :             {
    1208           0 :                 dfLatitude = 0.0;
    1209             :                 // For the most part, this is an invalid LCP, but it is a
    1210             :                 // changeable value in Flammap/Farsite, etc.  We should
    1211             :                 // probably be strict here all the time.
    1212           0 :                 CPLError(CE_Failure, CPLE_AppDefined,
    1213             :                          "Could not calculate latitude from spatial "
    1214             :                          "reference and LATITUDE was not set.");
    1215           0 :                 return nullptr;
    1216             :             }
    1217             :         }
    1218           0 :         OGRCoordinateTransformation::DestroyCT(poCT);
    1219             :     }
    1220             :     else
    1221             :     {
    1222             :         // See comment above on failure to transform.
    1223           0 :         CPLError(CE_Failure, CPLE_AppDefined,
    1224             :                  "Could not calculate latitude from spatial reference "
    1225             :                  "and LATITUDE was not set.");
    1226           0 :         return nullptr;
    1227             :     }
    1228             :     // Set the linear units if the metadata item was not already set, and we
    1229             :     // have an SRS.
    1230          33 :     if (bSetLinearUnits && poSrcSRS)
    1231             :     {
    1232           0 :         const char *pszUnit = poSrcSRS->GetAttrValue("UNIT", 0);
    1233           0 :         if (pszUnit == nullptr)
    1234             :         {
    1235           0 :             if (bStrict)
    1236             :             {
    1237           0 :                 CPLError(CE_Failure, CPLE_AppDefined,
    1238             :                          "Could not parse linear unit.");
    1239           0 :                 return nullptr;
    1240             :             }
    1241             :             else
    1242             :             {
    1243           0 :                 CPLError(CE_Warning, CPLE_AppDefined,
    1244             :                          "Could not parse linear unit, using meters");
    1245           0 :                 nLinearUnits = 0;
    1246             :             }
    1247             :         }
    1248             :         else
    1249             :         {
    1250           0 :             CPLDebug("LCP", "Setting linear unit to %s", pszUnit);
    1251           0 :             if (EQUAL(pszUnit, "meter") || EQUAL(pszUnit, "metre"))
    1252             :             {
    1253           0 :                 nLinearUnits = 0;
    1254             :             }
    1255           0 :             else if (EQUAL(pszUnit, "feet") || EQUAL(pszUnit, "foot"))
    1256             :             {
    1257           0 :                 nLinearUnits = 1;
    1258             :             }
    1259           0 :             else if (STARTS_WITH_CI(pszUnit, "kilomet"))
    1260             :             {
    1261           0 :                 nLinearUnits = 2;
    1262             :             }
    1263             :             else
    1264             :             {
    1265           0 :                 if (bStrict)
    1266           0 :                     nLinearUnits = 0;
    1267             :             }
    1268           0 :             pszUnit = poSrcSRS->GetAttrValue("UNIT", 1);
    1269           0 :             if (pszUnit != nullptr)
    1270             :             {
    1271           0 :                 const double dfScale = CPLAtof(pszUnit);
    1272           0 :                 if (dfScale != 1.0)
    1273             :                 {
    1274           0 :                     if (bStrict)
    1275             :                     {
    1276           0 :                         CPLError(CE_Failure, CPLE_AppDefined,
    1277             :                                  "Unit scale is %lf (!=1.0). It is not "
    1278             :                                  "supported.",
    1279             :                                  dfScale);
    1280           0 :                         return nullptr;
    1281             :                     }
    1282             :                     else
    1283             :                     {
    1284           0 :                         CPLError(CE_Warning, CPLE_AppDefined,
    1285             :                                  "Unit scale is %lf (!=1.0). It is not "
    1286             :                                  "supported, ignoring.",
    1287             :                                  dfScale);
    1288             :                     }
    1289             :                 }
    1290             :             }
    1291           0 :         }
    1292             :     }
    1293          33 :     else if (bSetLinearUnits)
    1294             :     {
    1295             :         // This can be defaulted if it is not a strict creation.
    1296           0 :         if (bStrict)
    1297             :         {
    1298           0 :             CPLError(CE_Failure, CPLE_AppDefined,
    1299             :                      "Could not parse linear unit from spatial reference "
    1300             :                      "and LINEAR_UNIT was not set.");
    1301           0 :             return nullptr;
    1302             :         }
    1303             :         else
    1304             :         {
    1305           0 :             CPLError(CE_Warning, CPLE_AppDefined,
    1306             :                      "Could not parse linear unit from spatial reference "
    1307             :                      "and LINEAR_UNIT was not set, defaulting to meters.");
    1308           0 :             nLinearUnits = 0;
    1309             :         }
    1310             :     }
    1311             : 
    1312          33 :     const char *pszDescription = CSLFetchNameValueDef(
    1313             :         papszOptions, "DESCRIPTION", "LCP file created by GDAL.");
    1314             : 
    1315             :     // Loop through and get the stats for the bands if we need to calculate
    1316             :     // them.  This probably should be done when we copy the data over to the
    1317             :     // destination dataset, since we load the values into memory, but this is
    1318             :     // much simpler code using GDALRasterBand->GetStatistics().  We also may
    1319             :     // need to classify the data (number of unique values and a list of those
    1320             :     // values if the number of unique values is > 100.  It is currently unclear
    1321             :     // how these data are used though, so we will implement that at some point
    1322             :     // if need be.
    1323          33 :     double *padfMin = static_cast<double *>(CPLMalloc(sizeof(double) * nBands));
    1324          33 :     double *padfMax = static_cast<double *>(CPLMalloc(sizeof(double) * nBands));
    1325             : 
    1326             :     // Initialize these arrays to zeros
    1327             :     GInt32 *panFound =
    1328          33 :         static_cast<GInt32 *>(VSIMalloc2(sizeof(GInt32), nBands));
    1329          33 :     memset(panFound, 0, sizeof(GInt32) * nBands);
    1330             :     GInt32 *panClasses = static_cast<GInt32 *>(
    1331          33 :         VSIMalloc3(sizeof(GInt32), nBands, LCP_MAX_CLASSES));
    1332          33 :     memset(panClasses, 0, sizeof(GInt32) * nBands * LCP_MAX_CLASSES);
    1333             : 
    1334          33 :     if (bCalculateStats)
    1335             :     {
    1336             : 
    1337         353 :         for (int i = 0; i < nBands; i++)
    1338             :         {
    1339         320 :             GDALRasterBand *poBand = poSrcDS->GetRasterBand(i + 1);
    1340         320 :             double dfDummy = 0.0;
    1341         640 :             CPLErr eErr = poBand->GetStatistics(
    1342         320 :                 FALSE, TRUE, &padfMin[i], &padfMax[i], &dfDummy, &dfDummy);
    1343         320 :             if (eErr != CE_None)
    1344             :             {
    1345           0 :                 CPLError(CE_Warning, CPLE_AppDefined,
    1346             :                          "Failed to properly "
    1347             :                          "calculate statistics "
    1348             :                          "on band %d",
    1349             :                          i);
    1350           0 :                 padfMin[i] = 0.0;
    1351           0 :                 padfMax[i] = 0.0;
    1352             :             }
    1353             : 
    1354             :             // See comment above.
    1355         320 :             if (bClassifyData)
    1356             :             {
    1357         640 :                 eErr = ClassifyBandData(poBand, panFound[i],
    1358         320 :                                         panClasses + (i * LCP_MAX_CLASSES));
    1359         320 :                 if (eErr != CE_None)
    1360             :                 {
    1361           0 :                     CPLError(CE_Warning, CPLE_AppDefined,
    1362             :                              "Failed to classify band data on band %d.", i);
    1363             :                 }
    1364             :             }
    1365             :         }
    1366             :     }
    1367             : 
    1368          33 :     VSILFILE *fp = VSIFOpenL(pszFilename, "wb");
    1369          33 :     if (fp == nullptr)
    1370             :     {
    1371           0 :         CPLError(CE_Failure, CPLE_OpenFailed, "Unable to create lcp file %s.",
    1372             :                  pszFilename);
    1373           0 :         CPLFree(padfMin);
    1374           0 :         CPLFree(padfMax);
    1375           0 :         CPLFree(panFound);
    1376           0 :         CPLFree(panClasses);
    1377           0 :         return nullptr;
    1378             :     }
    1379             : 
    1380             :     /* -------------------------------------------------------------------- */
    1381             :     /*      Write the header                                                */
    1382             :     /* -------------------------------------------------------------------- */
    1383             : 
    1384          33 :     GInt32 nTemp = bHaveCrownFuels ? 21 : 20;
    1385          33 :     CPL_LSBPTR32(&nTemp);
    1386          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&nTemp, 4, 1, fp));
    1387          33 :     nTemp = bHaveGroundFuels ? 21 : 20;
    1388          33 :     CPL_LSBPTR32(&nTemp);
    1389          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&nTemp, 4, 1, fp));
    1390             : 
    1391          33 :     const int nXSize = poSrcDS->GetRasterXSize();
    1392          33 :     nTemp = static_cast<GInt32>(dfLatitude + 0.5);
    1393          33 :     CPL_LSBPTR32(&nTemp);
    1394          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&nTemp, 4, 1, fp));
    1395          33 :     dfLongitude = srcGT[0] + srcGT[1] * nXSize;
    1396          33 :     CPL_LSBPTR64(&dfLongitude);
    1397          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&dfLongitude, 8, 1, fp));
    1398          33 :     dfLongitude = srcGT[0];
    1399          33 :     CPL_LSBPTR64(&dfLongitude);
    1400          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&dfLongitude, 8, 1, fp));
    1401          33 :     dfLatitude = srcGT[3];
    1402          33 :     CPL_LSBPTR64(&dfLatitude);
    1403          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&dfLatitude, 8, 1, fp));
    1404          33 :     dfLatitude = srcGT[3] + srcGT[5] * nYSize;
    1405          33 :     CPL_LSBPTR64(&dfLatitude);
    1406          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&dfLatitude, 8, 1, fp));
    1407             : 
    1408             :     // Swap the two classification arrays if we are writing them, and they need
    1409             :     // to be swapped.
    1410             : #ifdef CPL_MSB
    1411             :     if (bClassifyData)
    1412             :     {
    1413             :         GDALSwapWords(panFound, 2, nBands, 2);
    1414             :         GDALSwapWords(panClasses, 2, LCP_MAX_CLASSES, 2);
    1415             :     }
    1416             : #endif
    1417             : 
    1418          33 :     if (bCalculateStats)
    1419             :     {
    1420         353 :         for (int i = 0; i < nBands; i++)
    1421             :         {
    1422             :             // If we don't have Crown fuels, but do have Ground fuels, we
    1423             :             // have to 'fast forward'.
    1424         320 :             if (i == 5 && !bHaveCrownFuels && bHaveGroundFuels)
    1425             :             {
    1426           1 :                 CPL_IGNORE_RET_VAL(VSIFSeekL(fp, 3340, SEEK_SET));
    1427             :             }
    1428         320 :             nTemp = static_cast<GInt32>(padfMin[i]);
    1429         320 :             CPL_LSBPTR32(&nTemp);
    1430         320 :             CPL_IGNORE_RET_VAL(VSIFWriteL(&nTemp, 4, 1, fp));
    1431         320 :             nTemp = static_cast<GInt32>(padfMax[i]);
    1432         320 :             CPL_LSBPTR32(&nTemp);
    1433         320 :             CPL_IGNORE_RET_VAL(VSIFWriteL(&nTemp, 4, 1, fp));
    1434         320 :             if (bClassifyData)
    1435             :             {
    1436             :                 // These two arrays were swapped in their entirety above.
    1437         320 :                 CPL_IGNORE_RET_VAL(VSIFWriteL(panFound + i, 4, 1, fp));
    1438         320 :                 CPL_IGNORE_RET_VAL(
    1439         320 :                     VSIFWriteL(panClasses + (i * LCP_MAX_CLASSES), 4,
    1440             :                                LCP_MAX_CLASSES, fp));
    1441             :             }
    1442             :             else
    1443             :             {
    1444           0 :                 nTemp = -1;
    1445           0 :                 CPL_LSBPTR32(&nTemp);
    1446           0 :                 CPL_IGNORE_RET_VAL(VSIFWriteL(&nTemp, 4, 1, fp));
    1447           0 :                 CPL_IGNORE_RET_VAL(
    1448           0 :                     VSIFSeekL(fp, 4 * LCP_MAX_CLASSES, SEEK_CUR));
    1449             :             }
    1450             :         }
    1451             :     }
    1452             :     else
    1453             :     {
    1454           0 :         CPL_IGNORE_RET_VAL(VSIFSeekL(fp, 4164, SEEK_SET));
    1455             :     }
    1456          33 :     CPLFree(padfMin);
    1457          33 :     CPLFree(padfMax);
    1458          33 :     CPLFree(panFound);
    1459          33 :     CPLFree(panClasses);
    1460             : 
    1461             :     // Should be at one of 3 locations, 2104, 3340, or 4164.
    1462          33 :     CPLAssert(VSIFTellL(fp) == 2104 || VSIFTellL(fp) == 3340 ||
    1463             :               VSIFTellL(fp) == 4164);
    1464          33 :     CPL_IGNORE_RET_VAL(VSIFSeekL(fp, 4164, SEEK_SET));
    1465             : 
    1466             :     /* Image size */
    1467          33 :     nTemp = static_cast<GInt32>(nXSize);
    1468          33 :     CPL_LSBPTR32(&nTemp);
    1469          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&nTemp, 4, 1, fp));
    1470          33 :     nTemp = static_cast<GInt32>(nYSize);
    1471          33 :     CPL_LSBPTR32(&nTemp);
    1472          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&nTemp, 4, 1, fp));
    1473             : 
    1474             :     // X and Y boundaries.
    1475             :     // Max x.
    1476          33 :     double dfTemp = srcGT[0] + srcGT[1] * nXSize;
    1477          33 :     CPL_LSBPTR64(&dfTemp);
    1478          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&dfTemp, 8, 1, fp));
    1479             :     // Min x.
    1480          33 :     dfTemp = srcGT[0];
    1481          33 :     CPL_LSBPTR64(&dfTemp);
    1482          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&dfTemp, 8, 1, fp));
    1483             :     // Max y.
    1484          33 :     dfTemp = srcGT[3];
    1485          33 :     CPL_LSBPTR64(&dfTemp);
    1486          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&dfTemp, 8, 1, fp));
    1487             :     // Min y.
    1488          33 :     dfTemp = srcGT[3] + srcGT[5] * nYSize;
    1489          33 :     CPL_LSBPTR64(&dfTemp);
    1490          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&dfTemp, 8, 1, fp));
    1491             : 
    1492          33 :     nTemp = nLinearUnits;
    1493          33 :     CPL_LSBPTR32(&nTemp);
    1494          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&nTemp, 4, 1, fp));
    1495             : 
    1496             :     // Resolution.
    1497             :     // X resolution.
    1498          33 :     dfTemp = srcGT[1];
    1499          33 :     CPL_LSBPTR64(&dfTemp);
    1500          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&dfTemp, 8, 1, fp));
    1501             :     // Y resolution.
    1502          33 :     dfTemp = fabs(srcGT[5]);
    1503          33 :     CPL_LSBPTR64(&dfTemp);
    1504          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(&dfTemp, 8, 1, fp));
    1505             : 
    1506             : #ifdef CPL_MSB
    1507             :     GDALSwapWords(panMetadata, 2, LCP_MAX_BANDS, 2);
    1508             : #endif
    1509          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(panMetadata, 2, LCP_MAX_BANDS, fp));
    1510             : 
    1511             :     // Write the source filenames.
    1512          33 :     char **papszFileList = poSrcDS->GetFileList();
    1513          33 :     if (papszFileList != nullptr)
    1514             :     {
    1515         308 :         for (int i = 0; i < nBands; i++)
    1516             :         {
    1517         280 :             if (i == 5 && !bHaveCrownFuels && bHaveGroundFuels)
    1518             :             {
    1519           0 :                 CPL_IGNORE_RET_VAL(VSIFSeekL(fp, 6292, SEEK_SET));
    1520             :             }
    1521         280 :             CPL_IGNORE_RET_VAL(
    1522         280 :                 VSIFWriteL(papszFileList[0], 1,
    1523             :                            CPLStrnlen(papszFileList[0], LCP_MAX_PATH), fp));
    1524         280 :             CPL_IGNORE_RET_VAL(VSIFSeekL(fp, 4244 + (256 * (i + 1)), SEEK_SET));
    1525             :         }
    1526             :     }
    1527             :     // No file list, mem driver, etc.
    1528             :     else
    1529             :     {
    1530           5 :         CPL_IGNORE_RET_VAL(VSIFSeekL(fp, 6804, SEEK_SET));
    1531             :     }
    1532          33 :     CSLDestroy(papszFileList);
    1533             :     // Should be at location 5524, 6292 or 6804.
    1534          33 :     CPLAssert(VSIFTellL(fp) == 5524 || VSIFTellL(fp) == 6292 ||
    1535             :               VSIFTellL(fp) == 6804);
    1536          33 :     CPL_IGNORE_RET_VAL(VSIFSeekL(fp, 6804, SEEK_SET));
    1537             : 
    1538             :     // Description .
    1539          33 :     CPL_IGNORE_RET_VAL(VSIFWriteL(
    1540             :         pszDescription, 1, CPLStrnlen(pszDescription, LCP_MAX_DESC), fp));
    1541             : 
    1542             :     // Should be at or below location 7316, all done with the header.
    1543          33 :     CPLAssert(VSIFTellL(fp) <= 7316);
    1544          33 :     CPL_IGNORE_RET_VAL(VSIFSeekL(fp, 7316, SEEK_SET));
    1545             : 
    1546             :     /* -------------------------------------------------------------------- */
    1547             :     /*      Loop over image, copying image data.                            */
    1548             :     /* -------------------------------------------------------------------- */
    1549             : 
    1550          33 :     GInt16 *panScanline = static_cast<GInt16 *>(VSIMalloc3(2, nBands, nXSize));
    1551             : 
    1552          33 :     if (!pfnProgress(0.0, nullptr, pProgressData))
    1553             :     {
    1554           0 :         CPL_IGNORE_RET_VAL(VSIFCloseL(fp));
    1555           0 :         VSIFree(panScanline);
    1556           0 :         return nullptr;
    1557             :     }
    1558         399 :     for (int iLine = 0; iLine < nYSize; iLine++)
    1559             :     {
    1560        3826 :         for (int iBand = 0; iBand < nBands; iBand++)
    1561             :         {
    1562        3460 :             GDALRasterBand *poBand = poSrcDS->GetRasterBand(iBand + 1);
    1563        6920 :             CPLErr eErr = poBand->RasterIO(
    1564        3460 :                 GF_Read, 0, iLine, nXSize, 1, panScanline + iBand, nXSize, 1,
    1565        3460 :                 GDT_Int16, nBands * 2, static_cast<size_t>(nBands) * nXSize * 2,
    1566             :                 nullptr);
    1567             :             // Not sure what to do here.
    1568        3460 :             if (eErr != CE_None)
    1569             :             {
    1570           0 :                 CPLError(CE_Warning, CPLE_AppDefined,
    1571             :                          "Error reported in "
    1572             :                          "RasterIO");
    1573             :             }
    1574             :         }
    1575             : #ifdef CPL_MSB
    1576             :         GDALSwapWordsEx(panScanline, 2, static_cast<size_t>(nBands) * nXSize,
    1577             :                         2);
    1578             : #endif
    1579         366 :         CPL_IGNORE_RET_VAL(VSIFWriteL(
    1580         366 :             panScanline, 2, static_cast<size_t>(nBands) * nXSize, fp));
    1581             : 
    1582         366 :         if (!pfnProgress(iLine / static_cast<double>(nYSize), nullptr,
    1583             :                          pProgressData))
    1584             :         {
    1585           0 :             VSIFree(reinterpret_cast<void *>(panScanline));
    1586           0 :             CPL_IGNORE_RET_VAL(VSIFCloseL(fp));
    1587           0 :             return nullptr;
    1588             :         }
    1589             :     }
    1590          33 :     VSIFree(panScanline);
    1591          33 :     CPL_IGNORE_RET_VAL(VSIFCloseL(fp));
    1592          33 :     if (!pfnProgress(1.0, nullptr, pProgressData))
    1593             :     {
    1594           0 :         return nullptr;
    1595             :     }
    1596             : 
    1597             :     // Try to write projection file.  *Most* landfire data follows ESRI
    1598             :     // style projection files, so we use the same code as the AAIGrid driver.
    1599          33 :     if (poSrcSRS)
    1600             :     {
    1601           0 :         char *pszESRIProjection = nullptr;
    1602           0 :         const char *const apszOptions[] = {"FORMAT=WKT1_ESRI", nullptr};
    1603           0 :         poSrcSRS->exportToWkt(&pszESRIProjection, apszOptions);
    1604           0 :         if (pszESRIProjection)
    1605             :         {
    1606             :             char *const pszDirname =
    1607           0 :                 CPLStrdup(CPLGetPathSafe(pszFilename).c_str());
    1608             :             char *const pszBasename =
    1609           0 :                 CPLStrdup(CPLGetBasenameSafe(pszFilename).c_str());
    1610           0 :             char *pszPrjFilename = CPLStrdup(
    1611           0 :                 CPLFormFilenameSafe(pszDirname, pszBasename, "prj").c_str());
    1612           0 :             fp = VSIFOpenL(pszPrjFilename, "wt");
    1613           0 :             if (fp != nullptr)
    1614             :             {
    1615           0 :                 CPL_IGNORE_RET_VAL(VSIFWriteL(pszESRIProjection, 1,
    1616             :                                               strlen(pszESRIProjection), fp));
    1617           0 :                 CPL_IGNORE_RET_VAL(VSIFCloseL(fp));
    1618             :             }
    1619             :             else
    1620             :             {
    1621           0 :                 CPLError(CE_Failure, CPLE_FileIO, "Unable to create file %s.",
    1622             :                          pszPrjFilename);
    1623             :             }
    1624           0 :             CPLFree(pszDirname);
    1625           0 :             CPLFree(pszBasename);
    1626           0 :             CPLFree(pszPrjFilename);
    1627             :         }
    1628           0 :         CPLFree(pszESRIProjection);
    1629             :     }
    1630          33 :     return static_cast<GDALDataset *>(GDALOpen(pszFilename, GA_ReadOnly));
    1631             : }
    1632             : 
    1633             : /************************************************************************/
    1634             : /*                         GDALRegister_LCP()                           */
    1635             : /************************************************************************/
    1636             : 
    1637        2033 : void GDALRegister_LCP()
    1638             : 
    1639             : {
    1640        2033 :     if (GDALGetDriverByName("LCP") != nullptr)
    1641         283 :         return;
    1642             : 
    1643        1750 :     GDALDriver *poDriver = new GDALDriver();
    1644             : 
    1645        1750 :     poDriver->SetDescription("LCP");
    1646        1750 :     poDriver->SetMetadataItem(GDAL_DCAP_RASTER, "YES");
    1647        1750 :     poDriver->SetMetadataItem(GDAL_DMD_LONGNAME,
    1648        1750 :                               "FARSITE v.4 Landscape File (.lcp)");
    1649        1750 :     poDriver->SetMetadataItem(GDAL_DMD_EXTENSION, "lcp");
    1650        1750 :     poDriver->SetMetadataItem(GDAL_DMD_HELPTOPIC, "drivers/raster/lcp.html");
    1651             : 
    1652        1750 :     poDriver->SetMetadataItem(GDAL_DCAP_VIRTUALIO, "YES");
    1653             : 
    1654        1750 :     poDriver->SetMetadataItem(GDAL_DMD_CREATIONDATATYPES, "Int16");
    1655        1750 :     poDriver->SetMetadataItem(
    1656             :         GDAL_DMD_CREATIONOPTIONLIST,
    1657             :         "<CreationOptionList>"
    1658             :         "   <Option name='ELEVATION_UNIT' type='string-select' "
    1659             :         "default='METERS' description='Elevation units'>"
    1660             :         "       <Value>METERS</Value>"
    1661             :         "       <Value>FEET</Value>"
    1662             :         "   </Option>"
    1663             :         "   <Option name='SLOPE_UNIT' type='string-select' default='DEGREES' "
    1664             :         "description='Slope units'>"
    1665             :         "       <Value>DEGREES</Value>"
    1666             :         "       <Value>PERCENT</Value>"
    1667             :         "   </Option>"
    1668             :         "   <Option name='ASPECT_UNIT' type='string-select' "
    1669             :         "default='AZIMUTH_DEGREES'>"
    1670             :         "       <Value>GRASS_CATEGORIES</Value>"
    1671             :         "       <Value>AZIMUTH_DEGREES</Value>"
    1672             :         "       <Value>GRASS_DEGREES</Value>"
    1673             :         "   </Option>"
    1674             :         "   <Option name='FUEL_MODEL_OPTION' type='string-select' "
    1675             :         "default='NO_CUSTOM_AND_NO_FILE'>"
    1676             :         "       <Value>NO_CUSTOM_AND_NO_FILE</Value>"
    1677             :         "       <Value>CUSTOM_AND_NO_FILE</Value>"
    1678             :         "       <Value>NO_CUSTOM_AND_FILE</Value>"
    1679             :         "       <Value>CUSTOM_AND_FILE</Value>"
    1680             :         "   </Option>"
    1681             :         "   <Option name='CANOPY_COV_UNIT' type='string-select' "
    1682             :         "default='PERCENT'>"
    1683             :         "       <Value>CATEGORIES</Value>"
    1684             :         "       <Value>PERCENT</Value>"
    1685             :         "   </Option>"
    1686             :         "   <Option name='CANOPY_HT_UNIT' type='string-select' "
    1687             :         "default='METERS_X_10'>"
    1688             :         "       <Value>METERS</Value>"
    1689             :         "       <Value>FEET</Value>"
    1690             :         "       <Value>METERS_X_10</Value>"
    1691             :         "       <Value>FEET_X_10</Value>"
    1692             :         "   </Option>"
    1693             :         "   <Option name='CBH_UNIT' type='string-select' default='METERS_X_10'>"
    1694             :         "       <Value>METERS</Value>"
    1695             :         "       <Value>FEET</Value>"
    1696             :         "       <Value>METERS_X_10</Value>"
    1697             :         "       <Value>FEET_X_10</Value>"
    1698             :         "   </Option>"
    1699             :         "   <Option name='CBD_UNIT' type='string-select' "
    1700             :         "default='KG_PER_CUBIC_METER_X_100'>"
    1701             :         "       <Value>KG_PER_CUBIC_METER</Value>"
    1702             :         "       <Value>POUND_PER_CUBIC_FOOT</Value>"
    1703             :         "       <Value>KG_PER_CUBIC_METER_X_100</Value>"
    1704             :         "       <Value>POUND_PER_CUBIC_FOOT_X_1000</Value>"
    1705             :         "   </Option>"
    1706             :         "   <Option name='DUFF_UNIT' type='string-select' "
    1707             :         "default='MG_PER_HECTARE_X_10'>"
    1708             :         "       <Value>MG_PER_HECTARE_X_10</Value>"
    1709             :         "       <Value>TONS_PER_ACRE_X_10</Value>"
    1710             :         "   </Option>"
    1711             :         // Kyle does not think we need to override this, but maybe?
    1712             :         // "   <Option name='CWD_OPTION' type='boolean' default='FALSE'
    1713             :         // description='Override logic for setting the coarse woody presence'/>"
    1714             :         // */
    1715             :         "   <Option name='CALCULATE_STATS' type='boolean' default='YES' "
    1716             :         "description='Write the stats to the lcp'/>"
    1717             :         "   <Option name='CLASSIFY_DATA' type='boolean' default='YES' "
    1718             :         "description='Write the stats to the lcp'/>"
    1719             :         "   <Option name='LINEAR_UNIT' type='string-select' "
    1720             :         "default='SET_FROM_SRS' description='Set the linear units in the lcp'>"
    1721             :         "       <Value>SET_FROM_SRS</Value>"
    1722             :         "       <Value>METER</Value>"
    1723             :         "       <Value>FOOT</Value>"
    1724             :         "       <Value>KILOMETER</Value>"
    1725             :         "   </Option>"
    1726             :         "   <Option name='LATITUDE' type='int' default='' description='Set the "
    1727             :         "latitude for the dataset, this overrides the driver trying to set it "
    1728             :         "programmatically in EPSG:4269'/>"
    1729             :         "   <Option name='DESCRIPTION' type='string' default='LCP file created "
    1730             :         "by GDAL' description='A short description of the lcp file'/>"
    1731        1750 :         "</CreationOptionList>");
    1732             : 
    1733        1750 :     poDriver->pfnOpen = LCPDataset::Open;
    1734        1750 :     poDriver->pfnCreateCopy = LCPDataset::CreateCopy;
    1735        1750 :     poDriver->pfnIdentify = LCPDataset::Identify;
    1736             : 
    1737        1750 :     GetGDALDriverManager()->RegisterDriver(poDriver);
    1738             : }

Generated by: LCOV version 1.14