LCOV - code coverage report
Current view: top level - frmts/mrf - mrf_band.cpp (source / functions) Hit Total Coverage
Test: gdal_filtered.info Lines: 487 754 64.6 %
Date: 2024-05-04 12:52:34 Functions: 31 48 64.6 %

          Line data    Source code
       1             : /*
       2             :  * Copyright (c) 2002-2012, California Institute of Technology.
       3             :  * All rights reserved.  Based on Government Sponsored Research under contracts
       4             :  * NAS7-1407 and/or NAS7-03001.
       5             :  *
       6             :  * Redistribution and use in source and binary forms, with or without
       7             :  * modification, are permitted provided that the following conditions are met:
       8             :  *   1. Redistributions of source code must retain the above copyright notice,
       9             :  * this list of conditions and the following disclaimer.
      10             :  *   2. Redistributions in binary form must reproduce the above copyright
      11             :  * notice, this list of conditions and the following disclaimer in the
      12             :  * documentation and/or other materials provided with the distribution.
      13             :  *   3. Neither the name of the California Institute of Technology (Caltech),
      14             :  * its operating division the Jet Propulsion Laboratory (JPL), the National
      15             :  * Aeronautics and Space Administration (NASA), nor the names of its
      16             :  * contributors may be used to endorse or promote products derived from this
      17             :  * software without specific prior written permission.
      18             :  *
      19             :  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
      20             :  * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
      21             :  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
      22             :  * ARE DISCLAIMED. IN NO EVENT SHALL THE CALIFORNIA INSTITUTE OF TECHNOLOGY BE
      23             :  * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
      24             :  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
      25             :  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
      26             :  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
      27             :  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
      28             :  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
      29             :  * POSSIBILITY OF SUCH DAMAGE.
      30             :  *
      31             :  * Copyright 2014-2021 Esri
      32             :  *
      33             :  * Licensed under the Apache License, Version 2.0 (the "License");
      34             :  * you may not use this file except in compliance with the License.
      35             :  * You may obtain a copy of the License at
      36             :  *
      37             :  * http://www.apache.org/licenses/LICENSE-2.0
      38             :  *
      39             :  * Unless required by applicable law or agreed to in writing, software
      40             :  * distributed under the License is distributed on an "AS IS" BASIS,
      41             :  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
      42             :  * See the License for the specific language governing permissions and
      43             :  * limitations under the License.
      44             :  */
      45             : 
      46             : /******************************************************************************
      47             :  *
      48             :  * Project:  Meta Raster File Format Driver Implementation, RasterBand
      49             :  * Purpose:  Implementation of MRF band
      50             :  *
      51             :  * Author:   Lucian Plesea, Lucian.Plesea jpl.nasa.gov, lplesea esri.com
      52             :  *
      53             :  ****************************************************************************/
      54             : 
      55             : #include "marfa.h"
      56             : #include "gdal_priv.h"
      57             : #include "ogr_srs_api.h"
      58             : #include "ogr_spatialref.h"
      59             : 
      60             : #include <vector>
      61             : #include <cassert>
      62             : #include <zlib.h>
      63             : #if defined(ZSTD_SUPPORT)
      64             : #include <zstd.h>
      65             : #endif
      66             : 
      67             : using namespace std::chrono;
      68             : 
      69             : NAMESPACE_MRF_START
      70             : 
      71             : // packs a block of a given type, with a stride
      72             : // Count is the number of items that need to be copied
      73             : // These are separate to allow for optimization
      74             : 
      75             : template <typename T>
      76          24 : static void cpy_stride_in(void *dst, void *src, int c, int stride)
      77             : {
      78          24 :     T *s = reinterpret_cast<T *>(src);
      79          24 :     T *d = reinterpret_cast<T *>(dst);
      80             : 
      81     6291480 :     while (c--)
      82             :     {
      83     6291460 :         *d++ = *s;
      84     6291460 :         s += stride;
      85             :     }
      86          24 : }
      87             : 
      88             : template <typename T>
      89          33 : static void cpy_stride_out(void *dst, void *src, int c, int stride)
      90             : {
      91          33 :     T *s = reinterpret_cast<T *>(src);
      92          33 :     T *d = reinterpret_cast<T *>(dst);
      93             : 
      94     8650780 :     while (c--)
      95             :     {
      96     8650750 :         *d = *s++;
      97     8650750 :         d += stride;
      98             :     }
      99          33 : }
     100             : 
     101             : // Does every value in the buffer have the same value, using strict comparison
     102             : template <typename T>
     103        5100 : inline int isAllVal(const T *b, size_t bytecount, double ndv)
     104             : {
     105        5100 :     T val = static_cast<T>(ndv);
     106        5100 :     size_t count = bytecount / sizeof(T);
     107     1843723 :     for (; count; --count)
     108             :     {
     109     1843723 :         if (*(b++) != val)
     110             :         {
     111        5093 :             return FALSE;
     112             :         }
     113             :     }
     114           7 :     return TRUE;
     115             : }
     116             : 
     117             : // Dispatcher based on gdal data type
     118        5100 : static int isAllVal(GDALDataType gt, void *b, size_t bytecount, double ndv)
     119             : {
     120             :     // Test to see if it has data
     121        5100 :     int isempty = false;
     122             : 
     123             :     // A case branch in a temporary macro, conversion from gdal enum to type
     124             : #define TEST_T(GType, T)                                                       \
     125             :     case GType:                                                                \
     126             :         isempty = isAllVal(reinterpret_cast<T *>(b), bytecount, ndv);          \
     127             :         break
     128             : 
     129        5100 :     switch (gt)
     130             :     {
     131        4941 :         TEST_T(GDT_Byte, GByte);
     132          28 :         TEST_T(GDT_UInt16, GUInt16);
     133          27 :         TEST_T(GDT_Int16, GInt16);
     134          26 :         TEST_T(GDT_UInt32, GUInt32);
     135          26 :         TEST_T(GDT_Int32, GInt32);
     136          26 :         TEST_T(GDT_Float32, float);
     137          26 :         TEST_T(GDT_Float64, double);
     138           0 :         default:
     139           0 :             break;
     140             :     }
     141             : #undef TEST_T
     142             : 
     143        5100 :     return isempty;
     144             : }
     145             : 
     146             : // Swap bytes in place, unconditional
     147             : // cppcheck-suppress constParameterReference
     148           0 : static void swab_buff(buf_mgr &src, const ILImage &img)
     149             : {
     150             :     size_t i;
     151           0 :     switch (GDALGetDataTypeSize(img.dt))
     152             :     {
     153           0 :         case 16:
     154             :         {
     155           0 :             short int *b = (short int *)src.buffer;
     156           0 :             for (i = src.size / 2; i; b++, i--)
     157           0 :                 *b = swab16(*b);
     158           0 :             break;
     159             :         }
     160           0 :         case 32:
     161             :         {
     162           0 :             int *b = (int *)src.buffer;
     163           0 :             for (i = src.size / 4; i; b++, i--)
     164           0 :                 *b = swab32(*b);
     165           0 :             break;
     166             :         }
     167           0 :         case 64:
     168             :         {
     169           0 :             long long *b = (long long *)src.buffer;
     170           0 :             for (i = src.size / 8; i; b++, i--)
     171           0 :                 *b = swab64(*b);
     172           0 :             break;
     173             :         }
     174             :     }
     175           0 : }
     176             : 
     177             : // Similar to compress2() but with flags to control zlib features
     178             : // Returns true if it worked
     179           7 : static int ZPack(const buf_mgr &src, buf_mgr &dst, int flags)
     180             : {
     181             :     z_stream stream;
     182             :     int err;
     183             : 
     184           7 :     memset(&stream, 0, sizeof(stream));
     185           7 :     stream.next_in = (Bytef *)src.buffer;
     186           7 :     stream.avail_in = (uInt)src.size;
     187           7 :     stream.next_out = (Bytef *)dst.buffer;
     188           7 :     stream.avail_out = (uInt)dst.size;
     189             : 
     190           7 :     int level = flags & ZFLAG_LMASK;
     191           7 :     if (level > 9)
     192           0 :         level = 9;
     193           7 :     if (level < 1)
     194           0 :         level = 1;
     195           7 :     int wb = MAX_WBITS;
     196             :     // if gz flag is set, ignore raw request
     197           7 :     if (flags & ZFLAG_GZ)
     198           0 :         wb += 16;
     199           7 :     else if (flags & ZFLAG_RAW)
     200           0 :         wb = -wb;
     201           7 :     int memlevel = 8;  // Good compromise
     202           7 :     int strategy = (flags & ZFLAG_SMASK) >> 6;
     203           7 :     if (strategy > 4)
     204           0 :         strategy = 0;
     205             : 
     206           7 :     err = deflateInit2(&stream, level, Z_DEFLATED, wb, memlevel, strategy);
     207           7 :     if (err != Z_OK)
     208           0 :         return err;
     209             : 
     210           7 :     err = deflate(&stream, Z_FINISH);
     211           7 :     if (err != Z_STREAM_END)
     212             :     {
     213           0 :         deflateEnd(&stream);
     214           0 :         return false;
     215             :     }
     216           7 :     dst.size = stream.total_out;
     217           7 :     err = deflateEnd(&stream);
     218           7 :     return err == Z_OK;
     219             : }
     220             : 
     221             : // Similar to uncompress() from zlib, accepts the ZFLAG_RAW
     222             : // Return true if it worked
     223           7 : static int ZUnPack(const buf_mgr &src, buf_mgr &dst, int flags)
     224             : {
     225             : 
     226             :     z_stream stream;
     227             :     int err;
     228             : 
     229           7 :     memset(&stream, 0, sizeof(stream));
     230           7 :     stream.next_in = (Bytef *)src.buffer;
     231           7 :     stream.avail_in = (uInt)src.size;
     232           7 :     stream.next_out = (Bytef *)dst.buffer;
     233           7 :     stream.avail_out = (uInt)dst.size;
     234             : 
     235             :     // 32 means autodetec gzip or zlib header, negative 15 is for raw
     236           7 :     int wb = (ZFLAG_RAW & flags) ? -MAX_WBITS : 32 + MAX_WBITS;
     237           7 :     err = inflateInit2(&stream, wb);
     238           7 :     if (err != Z_OK)
     239           0 :         return false;
     240             : 
     241           7 :     err = inflate(&stream, Z_FINISH);
     242           7 :     if (err != Z_STREAM_END)
     243             :     {
     244           0 :         inflateEnd(&stream);
     245           0 :         return false;
     246             :     }
     247           7 :     dst.size = stream.total_out;
     248           7 :     err = inflateEnd(&stream);
     249           7 :     return err == Z_OK;
     250             : }
     251             : 
     252             : /*
     253             :  * Deflates a buffer, extrasize is the available size in the buffer past the
     254             :  * input If the output fits past the data, it uses that area, otherwise it uses
     255             :  * a temporary buffer and copies the data over the input on return, returning a
     256             :  * pointer to it. The output size is returned in src.size Returns nullptr when
     257             :  * compression failed
     258             :  */
     259           7 : static void *DeflateBlock(buf_mgr &src, size_t extrasize, int flags)
     260             : {
     261             :     // The one we might need to allocate
     262           7 :     void *dbuff = nullptr;
     263           7 :     buf_mgr dst = {src.buffer + src.size, extrasize};
     264             : 
     265             :     // Allocate a temp buffer if there is not sufficient space,
     266             :     // We need to have a bit more than half the buffer available
     267           7 :     if (extrasize < (src.size + 64))
     268             :     {
     269           7 :         dst.size = src.size + 64;
     270           7 :         dbuff = VSIMalloc(dst.size);
     271           7 :         dst.buffer = (char *)dbuff;
     272           7 :         if (!dst.buffer)
     273           0 :             return nullptr;
     274             :     }
     275             : 
     276           7 :     if (!ZPack(src, dst, flags))
     277             :     {
     278           0 :         CPLFree(dbuff);  // Safe to call with NULL
     279           0 :         return nullptr;
     280             :     }
     281           7 :     if (dst.size > src.size)
     282             :     {
     283           0 :         CPLError(CE_Failure, CPLE_AppDefined,
     284             :                  "DeflateBlock(): dst.size > src.size");
     285           0 :         CPLFree(dbuff);  // Safe to call with NULL
     286           0 :         return nullptr;
     287             :     }
     288             : 
     289             :     // source size is used to hold the output size
     290           7 :     src.size = dst.size;
     291             :     // If we didn't allocate a buffer, the receiver can use it already
     292           7 :     if (!dbuff)
     293           0 :         return dst.buffer;
     294             : 
     295             :     // If we allocated a buffer, we need to copy the data to the input buffer.
     296           7 :     memcpy(src.buffer, dbuff, src.size);
     297           7 :     CPLFree(dbuff);
     298           7 :     return src.buffer;
     299             : }
     300             : 
     301             : #if defined(ZSTD_SUPPORT)
     302             : 
     303          11 : CPL_NOSANITIZE_UNSIGNED_INT_OVERFLOW static void rankfilter(buf_mgr &src,
     304             :                                                             size_t factor)
     305             : {
     306             :     // Arange bytes by rank
     307          11 :     if (factor > 1)
     308             :     {
     309           7 :         std::vector<char> tempb(src.size);
     310           7 :         char *d = tempb.data();
     311          34 :         for (size_t j = 0; j < factor; j++)
     312     7077920 :             for (size_t i = j; i < src.size; i += factor)
     313     7077890 :                 *d++ = src.buffer[i];
     314           7 :         memcpy(src.buffer, tempb.data(), src.size);
     315             :     }
     316             :     // byte delta
     317          11 :     auto p = reinterpret_cast<GByte *>(src.buffer);
     318          11 :     auto guard = p + src.size;
     319          11 :     GByte b(0);
     320     8126480 :     while (p < guard)
     321             :     {
     322     8126460 :         GByte temp = *p;
     323     8126460 :         *p -= b;
     324     8126460 :         b = temp;
     325     8126460 :         p++;
     326             :     }
     327          11 : }
     328             : 
     329           9 : CPL_NOSANITIZE_UNSIGNED_INT_OVERFLOW static void derank(buf_mgr &src,
     330             :                                                         size_t factor)
     331             : {
     332             :     // undo delta
     333           9 :     auto p = reinterpret_cast<GByte *>(src.buffer);
     334           9 :     auto guard = p + src.size;
     335           9 :     GByte b(0);
     336     7602180 :     while (p < guard)
     337             :     {
     338     7602180 :         b += *p;
     339     7602180 :         *p = b;
     340     7602180 :         p++;
     341             :     }
     342           9 :     if (factor > 1)
     343             :     {  // undo rank separation
     344           7 :         std::vector<char> tempb(src.size);
     345           7 :         char *d = tempb.data();
     346           7 :         size_t chunk = src.size / factor;
     347     1835020 :         for (size_t i = 0; i < chunk; i++)
     348     8912900 :             for (size_t j = 0; j < factor; j++)
     349     7077890 :                 *d++ = src.buffer[chunk * j + i];
     350           7 :         memcpy(src.buffer, tempb.data(), src.size);
     351             :     }
     352           9 : }
     353             : 
     354             : /*
     355             :  * Compress a buffer using zstd, extrasize is the available size in the buffer
     356             :  * past the input If ranks > 0, apply the rank filter If the output fits past
     357             :  * the data, it uses that area, otherwise it uses a temporary buffer and copies
     358             :  * the data over the input on return, returning a pointer to it. The output size
     359             :  * is returned in src.size Returns nullptr when compression failed
     360             :  */
     361          11 : static void *ZstdCompBlock(buf_mgr &src, size_t extrasize, int c_level,
     362             :                            ZSTD_CCtx *cctx, size_t ranks)
     363             : {
     364          11 :     if (!cctx)
     365           0 :         return nullptr;
     366          11 :     if (ranks && (src.size % ranks) == 0)
     367          11 :         rankfilter(src, ranks);
     368             : 
     369             :     // might need a buffer for the zstd output
     370          22 :     std::vector<char> dbuff;
     371          11 :     void *dst = src.buffer + src.size;
     372          11 :     size_t size = extrasize;
     373             :     // Allocate a temp buffer if there is not sufficient space.
     374             :     // Zstd bound is about (size * 1.004 + 64)
     375          11 :     if (size < ZSTD_compressBound(src.size))
     376             :     {
     377          11 :         size = ZSTD_compressBound(src.size);
     378          11 :         dbuff.resize(size);
     379          11 :         dst = dbuff.data();
     380             :     }
     381             : 
     382             :     // Use the streaming interface, it's faster and better
     383             :     // See discussion at https://github.com/facebook/zstd/issues/3729
     384          11 :     ZSTD_outBuffer output = {dst, size, 0};
     385          11 :     ZSTD_inBuffer input = {src.buffer, src.size, 0};
     386             :     // Set level
     387          11 :     ZSTD_CCtx_setParameter(cctx, ZSTD_c_compressionLevel, c_level);
     388             :     // First, pass a continue flag, otherwise it will compress in one go
     389          11 :     size_t val = ZSTD_compressStream2(cctx, &output, &input, ZSTD_e_continue);
     390             :     // If it worked, pass the end flag to flush the buffer
     391          11 :     if (val == 0)
     392          11 :         val = ZSTD_compressStream2(cctx, &output, &input, ZSTD_e_end);
     393          11 :     if (ZSTD_isError(val))
     394           0 :         return nullptr;
     395          11 :     val = output.pos;
     396             : 
     397             :     // If we didn't need the buffer, packed data is already in the user buffer
     398          11 :     if (dbuff.empty())
     399             :     {
     400           0 :         src.size = val;
     401           0 :         return dst;
     402             :     }
     403             : 
     404          11 :     if (val > (src.size + extrasize))
     405             :     {  // Doesn't fit in user buffer
     406           0 :         CPLError(CE_Failure, CPLE_AssertionFailed,
     407             :                  "MRF: ZSTD compression buffer too small");
     408           0 :         return nullptr;  // Error
     409             :     }
     410             : 
     411          11 :     memcpy(src.buffer, dbuff.data(), val);
     412          11 :     src.size = val;
     413          11 :     return src.buffer;
     414             : }
     415             : #endif
     416             : 
     417             : //
     418             : // The deflate_flags are available in all bands even if the DEFLATE option
     419             : // itself is not set.  This allows for PNG features to be controlled, as well
     420             : // as any other bands that use zlib by itself
     421             : //
     422         530 : MRFRasterBand::MRFRasterBand(MRFDataset *parent_dataset, const ILImage &image,
     423         530 :                              int band, int ov)
     424             :     : poMRFDS(parent_dataset),
     425         530 :       dodeflate(GetOptlist().FetchBoolean("DEFLATE", FALSE)),
     426             :       // Bring the quality to 0 to 9
     427         530 :       deflate_flags(image.quality / 10),
     428         530 :       dozstd(GetOptlist().FetchBoolean("ZSTD", FALSE)), zstd_level(9), m_l(ov),
     429        1590 :       img(image)
     430             : {
     431         530 :     nBand = band;
     432         530 :     eDataType = parent_dataset->current.dt;
     433         530 :     nRasterXSize = img.size.x;
     434         530 :     nRasterYSize = img.size.y;
     435         530 :     nBlockXSize = img.pagesize.x;
     436         530 :     nBlockYSize = img.pagesize.y;
     437         530 :     nBlocksPerRow = img.pagecount.x;
     438         530 :     nBlocksPerColumn = img.pagecount.y;
     439         530 :     img.NoDataValue = MRFRasterBand::GetNoDataValue(&img.hasNoData);
     440             : 
     441             :     // Pick up the twists, aka GZ, RAWZ headers
     442         530 :     if (GetOptlist().FetchBoolean("GZ", FALSE))
     443           0 :         deflate_flags |= ZFLAG_GZ;
     444         530 :     else if (GetOptlist().FetchBoolean("RAWZ", FALSE))
     445           0 :         deflate_flags |= ZFLAG_RAW;
     446             :     // And Pick up the ZLIB strategy, if any
     447         530 :     const char *zstrategy = GetOptlist().FetchNameValueDef("Z_STRATEGY", "");
     448         530 :     int zv = Z_DEFAULT_STRATEGY;
     449         530 :     if (EQUAL(zstrategy, "Z_HUFFMAN_ONLY"))
     450           0 :         zv = Z_HUFFMAN_ONLY;
     451         530 :     else if (EQUAL(zstrategy, "Z_RLE"))
     452           0 :         zv = Z_RLE;
     453         530 :     else if (EQUAL(zstrategy, "Z_FILTERED"))
     454           0 :         zv = Z_FILTERED;
     455         530 :     else if (EQUAL(zstrategy, "Z_FIXED"))
     456           0 :         zv = Z_FIXED;
     457         530 :     deflate_flags |= (zv << 6);
     458         530 :     if (image.quality < 23 && image.quality > 0)
     459           0 :         zstd_level = image.quality;
     460             : 
     461             : #if !defined(ZSTD_SUPPORT)
     462             :     if (dozstd)
     463             :     {  // signal error condition to caller
     464             :         CPLError(CE_Failure, CPLE_AssertionFailed,
     465             :                  "MRF: ZSTD support is not available");
     466             :         dozstd = FALSE;
     467             :     }
     468             : #endif
     469             :     // Chose zstd over deflate if both are enabled and available
     470         530 :     if (dozstd && dodeflate)
     471           0 :         dodeflate = FALSE;
     472         530 : }
     473             : 
     474             : // Clean up the overviews if they exist
     475        1060 : MRFRasterBand::~MRFRasterBand()
     476             : {
     477         602 :     while (!overviews.empty())
     478             :     {
     479          72 :         delete overviews.back();
     480          72 :         overviews.pop_back();
     481             :     }
     482         530 : }
     483             : 
     484             : // Look for a string from the dataset options or from the environment
     485          75 : const char *MRFRasterBand::GetOptionValue(const char *opt,
     486             :                                           const char *def) const
     487             : {
     488          75 :     const char *optValue = poMRFDS->optlist.FetchNameValue(opt);
     489          75 :     if (optValue)
     490           8 :         return optValue;
     491          67 :     return CPLGetConfigOption(opt, def);
     492             : }
     493             : 
     494             : // Utility function, returns a value from a vector corresponding to the band
     495             : // index or the first entry
     496         157 : static double getBandValue(const std::vector<double> &v, int idx)
     497             : {
     498         157 :     return (static_cast<int>(v.size()) > idx) ? v[idx] : v[0];
     499             : }
     500             : 
     501             : // Maybe we should check against the type range?
     502             : // It is not keeping track of how many values have been set,
     503             : // so the application should set none or all the bands
     504             : // This call is only valid during Create
     505          16 : CPLErr MRFRasterBand::SetNoDataValue(double val)
     506             : {
     507          16 :     if (poMRFDS->bCrystalized)
     508             :     {
     509           0 :         CPLError(CE_Failure, CPLE_AssertionFailed,
     510             :                  "MRF: NoData can be set only during file create");
     511           0 :         return CE_Failure;
     512             :     }
     513          16 :     if (GInt32(poMRFDS->vNoData.size()) < nBand)
     514           0 :         poMRFDS->vNoData.resize(nBand);
     515          16 :     poMRFDS->vNoData[nBand - 1] = val;
     516             :     // We also need to set it for this band
     517          16 :     img.NoDataValue = val;
     518          16 :     img.hasNoData = true;
     519          16 :     return CE_None;
     520             : }
     521             : 
     522        5973 : double MRFRasterBand::GetNoDataValue(int *pbSuccess)
     523             : {
     524        5973 :     const std::vector<double> &v = poMRFDS->vNoData;
     525        5973 :     if (v.empty())
     526        5816 :         return GDALPamRasterBand::GetNoDataValue(pbSuccess);
     527         157 :     if (pbSuccess)
     528         155 :         *pbSuccess = TRUE;
     529         157 :     return getBandValue(v, nBand - 1);
     530             : }
     531             : 
     532           0 : double MRFRasterBand::GetMinimum(int *pbSuccess)
     533             : {
     534           0 :     const std::vector<double> &v = poMRFDS->vMin;
     535           0 :     if (v.empty())
     536           0 :         return GDALPamRasterBand::GetMinimum(pbSuccess);
     537           0 :     if (pbSuccess)
     538           0 :         *pbSuccess = TRUE;
     539           0 :     return getBandValue(v, nBand - 1);
     540             : }
     541             : 
     542           0 : double MRFRasterBand::GetMaximum(int *pbSuccess)
     543             : {
     544           0 :     const std::vector<double> &v = poMRFDS->vMax;
     545           0 :     if (v.empty())
     546           0 :         return GDALPamRasterBand::GetMaximum(pbSuccess);
     547           0 :     if (pbSuccess)
     548           0 :         *pbSuccess = TRUE;
     549           0 :     return getBandValue(v, nBand - 1);
     550             : }
     551             : 
     552             : // Fill with typed ndv, count is always in bytes
     553             : template <typename T>
     554           0 : static CPLErr buff_fill(void *b, size_t count, const T ndv)
     555             : {
     556           0 :     T *buffer = static_cast<T *>(b);
     557           0 :     count /= sizeof(T);
     558           0 :     while (count--)
     559           0 :         *buffer++ = ndv;
     560           0 :     return CE_None;
     561             : }
     562             : 
     563             : /**
     564             :  *\brief Fills a buffer with no data
     565             :  *
     566             :  */
     567          56 : CPLErr MRFRasterBand::FillBlock(void *buffer)
     568             : {
     569             :     int success;
     570          56 :     double ndv = GetNoDataValue(&success);
     571          56 :     if (!success)
     572          56 :         ndv = 0.0;
     573          56 :     size_t bsb = blockSizeBytes();
     574             : 
     575             :     // use memset for speed for bytes, or if nodata is zeros
     576          56 :     if (eDataType == GDT_Byte || 0.0L == ndv)
     577             :     {
     578          56 :         memset(buffer, int(ndv), bsb);
     579          56 :         return CE_None;
     580             :     }
     581             : 
     582             : #define bf(T) buff_fill<T>(buffer, bsb, T(ndv));
     583           0 :     switch (eDataType)
     584             :     {
     585           0 :         case GDT_UInt16:
     586           0 :             return bf(GUInt16);
     587           0 :         case GDT_Int16:
     588           0 :             return bf(GInt16);
     589           0 :         case GDT_UInt32:
     590           0 :             return bf(GUInt32);
     591           0 :         case GDT_Int32:
     592           0 :             return bf(GInt32);
     593           0 :         case GDT_Float32:
     594           0 :             return bf(float);
     595           0 :         case GDT_Float64:
     596           0 :             return bf(double);
     597           0 :         default:
     598           0 :             break;
     599             :     }
     600             : #undef bf
     601             :     // Should exit before
     602           0 :     return CE_Failure;
     603             : }
     604             : 
     605             : /*\brief Interleave block fill
     606             :  *
     607             :  *  Acquire space for all the other bands, fill each one then drop the locks
     608             :  *  The current band output goes directly into the buffer
     609             :  */
     610             : 
     611           0 : CPLErr MRFRasterBand::FillBlock(int xblk, int yblk, void *buffer)
     612             : {
     613           0 :     std::vector<GDALRasterBlock *> blocks;
     614             : 
     615           0 :     for (int i = 0; i < poMRFDS->nBands; i++)
     616             :     {
     617           0 :         GDALRasterBand *b = poMRFDS->GetRasterBand(i + 1);
     618           0 :         if (b->GetOverviewCount() && 0 != m_l)
     619           0 :             b = b->GetOverview(m_l - 1);
     620             : 
     621             :         // Get the other band blocks, keep them around until later
     622           0 :         if (b == this)
     623             :         {
     624           0 :             FillBlock(buffer);
     625             :         }
     626             :         else
     627             :         {
     628           0 :             GDALRasterBlock *poBlock = b->GetLockedBlockRef(xblk, yblk, 1);
     629           0 :             if (poBlock == nullptr)  // Didn't get this block
     630           0 :                 break;
     631           0 :             FillBlock(poBlock->GetDataRef());
     632           0 :             blocks.push_back(poBlock);
     633             :         }
     634             :     }
     635             : 
     636             :     // Drop the locks for blocks we acquired
     637           0 :     for (int i = 0; i < int(blocks.size()); i++)
     638           0 :         blocks[i]->DropLock();
     639             : 
     640           0 :     return CE_None;
     641             : }
     642             : 
     643             : /*\brief Interleave block read
     644             :  *
     645             :  *  Acquire space for all the other bands, unpack from the dataset buffer, then
     646             :  * drop the locks The current band output goes directly into the buffer
     647             :  */
     648             : 
     649           8 : CPLErr MRFRasterBand::ReadInterleavedBlock(int xblk, int yblk, void *buffer)
     650             : {
     651           8 :     std::vector<GDALRasterBlock *> blocks;
     652             : 
     653          32 :     for (int i = 0; i < poMRFDS->nBands; i++)
     654             :     {
     655          24 :         GDALRasterBand *b = poMRFDS->GetRasterBand(i + 1);
     656          24 :         if (b->GetOverviewCount() && 0 != m_l)
     657           0 :             b = b->GetOverview(m_l - 1);
     658             : 
     659          24 :         void *ob = buffer;
     660             :         // Get the other band blocks, keep them around until later
     661          24 :         if (b != this)
     662             :         {
     663          16 :             GDALRasterBlock *poBlock = b->GetLockedBlockRef(xblk, yblk, 1);
     664          16 :             if (poBlock == nullptr)
     665           0 :                 break;
     666          16 :             ob = poBlock->GetDataRef();
     667          16 :             blocks.push_back(poBlock);
     668             :         }
     669             : 
     670             :         // Just the right mix of templates and macros make deinterleaving tidy
     671          24 :         void *pbuffer = poMRFDS->GetPBuffer();
     672             : #define CpySI(T)                                                               \
     673             :     cpy_stride_in<T>(ob, reinterpret_cast<T *>(pbuffer) + i,                   \
     674             :                      blockSizeBytes() / sizeof(T), img.pagesize.c)
     675             : 
     676             :         // Page is already in poMRFDS->pbuffer, not empty
     677             :         // There are only four cases, since only the data size matters
     678          24 :         switch (GDALGetDataTypeSize(eDataType) / 8)
     679             :         {
     680          24 :             case 1:
     681          24 :                 CpySI(GByte);
     682          24 :                 break;
     683           0 :             case 2:
     684           0 :                 CpySI(GInt16);
     685           0 :                 break;
     686           0 :             case 4:
     687           0 :                 CpySI(GInt32);
     688           0 :                 break;
     689           0 :             case 8:
     690           0 :                 CpySI(GIntBig);
     691           0 :                 break;
     692             :         }
     693             :     }
     694             : 
     695             : #undef CpySI
     696             : 
     697             :     // Drop the locks we acquired
     698          24 :     for (int i = 0; i < int(blocks.size()); i++)
     699          16 :         blocks[i]->DropLock();
     700             : 
     701          16 :     return CE_None;
     702             : }
     703             : 
     704             : /**
     705             :  *\brief Fetch a block from the backing store dataset and keep a copy in the
     706             :  *cache
     707             :  *
     708             :  * @param xblk The X block number, zero based
     709             :  * @param yblk The Y block number, zero based
     710             :  * @param buffer buffer
     711             :  *
     712             :  */
     713           4 : CPLErr MRFRasterBand::FetchBlock(int xblk, int yblk, void *buffer)
     714             : {
     715           4 :     assert(!poMRFDS->source.empty());
     716           4 :     CPLDebug("MRF_IB", "FetchBlock %d,%d,0,%d, level  %d\n", xblk, yblk, nBand,
     717             :              m_l);
     718             : 
     719           4 :     if (poMRFDS->clonedSource)  // This is a clone
     720           1 :         return FetchClonedBlock(xblk, yblk, buffer);
     721             : 
     722           3 :     const GInt32 cstride = img.pagesize.c;  // 1 if band separate
     723           3 :     ILSize req(xblk, yblk, 0, (nBand - 1) / cstride, m_l);
     724           3 :     GUIntBig infooffset = IdxOffset(req, img);
     725             : 
     726           3 :     GDALDataset *poSrcDS = nullptr;
     727           3 :     if (nullptr == (poSrcDS = poMRFDS->GetSrcDS()))
     728             :     {
     729           1 :         CPLError(CE_Failure, CPLE_AppDefined, "MRF: Can't open source file %s",
     730           1 :                  poMRFDS->source.c_str());
     731           1 :         return CE_Failure;
     732             :     }
     733             : 
     734             :     // Scale to base resolution
     735           2 :     double scl = pow(poMRFDS->scale, m_l);
     736           2 :     if (0 == m_l)
     737           2 :         scl = 1;  // To allow for precision issues
     738             : 
     739             :     // Prepare parameters for RasterIO, they might be different from a full page
     740           2 :     const GSpacing vsz = GDALGetDataTypeSizeBytes(eDataType);
     741           2 :     int Xoff = int(xblk * img.pagesize.x * scl + 0.5);
     742           2 :     int Yoff = int(yblk * img.pagesize.y * scl + 0.5);
     743           2 :     int readszx = int(img.pagesize.x * scl + 0.5);
     744           2 :     int readszy = int(img.pagesize.y * scl + 0.5);
     745             : 
     746             :     // Compare with the full size and clip to the right and bottom if needed
     747           2 :     int clip = 0;
     748           2 :     if (Xoff + readszx > poMRFDS->full.size.x)
     749             :     {
     750           2 :         clip |= 1;
     751           2 :         readszx = poMRFDS->full.size.x - Xoff;
     752             :     }
     753           2 :     if (Yoff + readszy > poMRFDS->full.size.y)
     754             :     {
     755           2 :         clip |= 1;
     756           2 :         readszy = poMRFDS->full.size.y - Yoff;
     757             :     }
     758             : 
     759             :     // This is where the whole page fits
     760           2 :     void *ob = buffer;
     761           2 :     if (cstride != 1)
     762           0 :         ob = poMRFDS->GetPBuffer();
     763             : 
     764             :     // Fill buffer with NoData if clipping
     765           2 :     if (clip)
     766           2 :         FillBlock(ob);
     767             : 
     768             :     // Use the dataset RasterIO to read one or all bands if interleaved
     769           4 :     CPLErr ret = poSrcDS->RasterIO(
     770             :         GF_Read, Xoff, Yoff, readszx, readszy, ob, pcount(readszx, int(scl)),
     771             :         pcount(readszy, int(scl)), eDataType, cstride,
     772           2 :         (1 == cstride) ? &nBand : nullptr, vsz * cstride,
     773           2 :         vsz * cstride * img.pagesize.x,
     774           2 :         (cstride != 1) ? vsz : vsz * img.pagesize.x * img.pagesize.y, nullptr);
     775             : 
     776           2 :     if (ret != CE_None)
     777           0 :         return ret;
     778             : 
     779             :     // Might have the block in the pbuffer, mark it anyhow
     780           2 :     poMRFDS->tile = req;
     781             :     buf_mgr filesrc;
     782           2 :     filesrc.buffer = (char *)ob;
     783           2 :     filesrc.size = static_cast<size_t>(img.pageSizeBytes);
     784             : 
     785           2 :     if (poMRFDS->bypass_cache)
     786             :     {  // No local caching, just return the data
     787           0 :         if (1 == cstride)
     788           0 :             return CE_None;
     789           0 :         return ReadInterleavedBlock(xblk, yblk, buffer);
     790             :     }
     791             : 
     792             :     // Test to see if it needs to be written, or just marked as checked
     793             :     int success;
     794           2 :     double val = GetNoDataValue(&success);
     795           2 :     if (!success)
     796           2 :         val = 0.0;
     797             : 
     798             :     // TODO: test band by band if data is interleaved
     799           2 :     if (isAllVal(eDataType, ob, img.pageSizeBytes, val))
     800             :     {
     801             :         // Mark it empty and checked, ignore the possible write error
     802           0 :         poMRFDS->WriteTile((void *)1, infooffset, 0);
     803           0 :         if (1 == cstride)
     804           0 :             return CE_None;
     805           0 :         return ReadInterleavedBlock(xblk, yblk, buffer);
     806             :     }
     807             : 
     808             :     // Write the page in the local cache
     809             : 
     810             :     // Have to use a separate buffer for compression output.
     811           2 :     void *outbuff = VSIMalloc(poMRFDS->pbsize);
     812           2 :     if (nullptr == outbuff)
     813             :     {
     814           0 :         CPLError(CE_Failure, CPLE_AppDefined,
     815             :                  "Can't get buffer for writing page");
     816             :         // This is not really an error for a cache, the data is fine
     817           0 :         return CE_Failure;
     818             :     }
     819             : 
     820           2 :     buf_mgr filedst = {static_cast<char *>(outbuff), poMRFDS->pbsize};
     821           2 :     auto start_time = steady_clock::now();
     822           2 :     Compress(filedst, filesrc);
     823             : 
     824             :     // Where the output is, in case we deflate
     825           2 :     void *usebuff = outbuff;
     826           2 :     if (dodeflate)
     827             :     {
     828           0 :         usebuff = DeflateBlock(filedst, poMRFDS->pbsize - filedst.size,
     829             :                                deflate_flags);
     830           0 :         if (!usebuff)
     831             :         {
     832           0 :             CPLError(CE_Failure, CPLE_AppDefined, "MRF: Deflate error");
     833           0 :             return CE_Failure;
     834             :         }
     835             :     }
     836             : 
     837             : #if defined(ZSTD_SUPPORT)
     838           2 :     else if (dozstd)
     839             :     {
     840           0 :         size_t ranks = 0;  // Assume no need for byte rank sort
     841           0 :         if (img.comp == IL_NONE || img.comp == IL_ZSTD)
     842           0 :             ranks =
     843           0 :                 static_cast<size_t>(GDALGetDataTypeSizeBytes(img.dt)) * cstride;
     844           0 :         usebuff = ZstdCompBlock(filedst, poMRFDS->pbsize - filedst.size,
     845           0 :                                 zstd_level, poMRFDS->getzsc(), ranks);
     846           0 :         if (!usebuff)
     847             :         {
     848           0 :             CPLError(CE_Failure, CPLE_AppDefined,
     849             :                      "MRF: ZSTD compression error");
     850           0 :             return CE_Failure;
     851             :         }
     852             :     }
     853             : #endif
     854             : 
     855           2 :     poMRFDS->write_timer +=
     856           2 :         duration_cast<nanoseconds>(steady_clock::now() - start_time);
     857             : 
     858             :     // Write and update the tile index
     859           2 :     ret = poMRFDS->WriteTile(usebuff, infooffset, filedst.size);
     860           2 :     CPLFree(outbuff);
     861             : 
     862             :     // If we hit an error or if unpaking is not needed
     863           2 :     if (ret != CE_None || cstride == 1)
     864           2 :         return ret;
     865             : 
     866             :     // data is already in DS buffer, deinterlace it in pixel blocks
     867           0 :     return ReadInterleavedBlock(xblk, yblk, buffer);
     868             : }
     869             : 
     870             : /**
     871             :  *\brief Fetch for a cloned MRF
     872             :  *
     873             :  * @param xblk The X block number, zero based
     874             :  * @param yblk The Y block number, zero based
     875             :  * @param buffer buffer
     876             :  *
     877             :  */
     878             : 
     879           1 : CPLErr MRFRasterBand::FetchClonedBlock(int xblk, int yblk, void *buffer)
     880             : {
     881           1 :     CPLDebug("MRF_IB", "FetchClonedBlock %d,%d,0,%d, level  %d\n", xblk, yblk,
     882             :              nBand, m_l);
     883             : 
     884             :     // Paranoid check
     885           1 :     assert(poMRFDS->clonedSource);
     886           1 :     MRFDataset *poSrc = static_cast<MRFDataset *>(poMRFDS->GetSrcDS());
     887           1 :     if (nullptr == poSrc)
     888             :     {
     889           0 :         CPLError(CE_Failure, CPLE_AppDefined, "MRF: Can't open source file %s",
     890           0 :                  poMRFDS->source.c_str());
     891           0 :         return CE_Failure;
     892             :     }
     893             : 
     894           1 :     if (poMRFDS->bypass_cache || GF_Read == DataMode())
     895             :     {
     896             :         // Can't store, so just fetch from source, which is an MRF with
     897             :         // identical structure
     898             :         MRFRasterBand *b =
     899           0 :             static_cast<MRFRasterBand *>(poSrc->GetRasterBand(nBand));
     900           0 :         if (b->GetOverviewCount() && m_l)
     901           0 :             b = static_cast<MRFRasterBand *>(b->GetOverview(m_l - 1));
     902           0 :         if (b == nullptr)
     903           0 :             return CE_Failure;
     904           0 :         return b->IReadBlock(xblk, yblk, buffer);
     905             :     }
     906             : 
     907           1 :     ILSize req(xblk, yblk, 0, (nBand - 1) / img.pagesize.c, m_l);
     908             :     ILIdx tinfo;
     909             : 
     910             :     // Get the cloned source tile info
     911             :     // The cloned source index is after the current one
     912           1 :     if (CE_None != poMRFDS->ReadTileIdx(tinfo, req, img, poMRFDS->idxSize))
     913             :     {
     914           0 :         CPLError(CE_Failure, CPLE_AppDefined,
     915             :                  "MRF: Unable to read cloned index entry");
     916           0 :         return CE_Failure;
     917             :     }
     918             : 
     919           1 :     GUIntBig infooffset = IdxOffset(req, img);
     920             :     CPLErr err;
     921             : 
     922             :     // Does the source have this tile?
     923           1 :     if (tinfo.size == 0)
     924             :     {  // Nope, mark it empty and return fill
     925           0 :         err = poMRFDS->WriteTile((void *)1, infooffset, 0);
     926           0 :         if (CE_None != err)
     927           0 :             return err;
     928           0 :         return FillBlock(buffer);
     929             :     }
     930             : 
     931           1 :     VSILFILE *srcfd = poSrc->DataFP();
     932           1 :     if (nullptr == srcfd)
     933             :     {
     934           0 :         CPLError(CE_Failure, CPLE_AppDefined,
     935             :                  "MRF: Can't open source data file %s",
     936           0 :                  poMRFDS->source.c_str());
     937           0 :         return CE_Failure;
     938             :     }
     939             : 
     940             :     // Need to read the tile from the source
     941           1 :     if (tinfo.size <= 0 || tinfo.size > INT_MAX)
     942             :     {
     943           0 :         CPLError(CE_Failure, CPLE_OutOfMemory,
     944             :                  "Invalid tile size " CPL_FRMT_GIB, tinfo.size);
     945           0 :         return CE_Failure;
     946             :     }
     947           1 :     char *buf = static_cast<char *>(VSIMalloc(static_cast<size_t>(tinfo.size)));
     948           1 :     if (buf == nullptr)
     949             :     {
     950           0 :         CPLError(CE_Failure, CPLE_OutOfMemory,
     951             :                  "Cannot allocate " CPL_FRMT_GIB " bytes", tinfo.size);
     952           0 :         return CE_Failure;
     953             :     }
     954             : 
     955           1 :     VSIFSeekL(srcfd, tinfo.offset, SEEK_SET);
     956           2 :     if (tinfo.size !=
     957           1 :         GIntBig(VSIFReadL(buf, 1, static_cast<size_t>(tinfo.size), srcfd)))
     958             :     {
     959           0 :         CPLFree(buf);
     960           0 :         CPLError(CE_Failure, CPLE_AppDefined,
     961             :                  "MRF: Can't read data from source %s",
     962             :                  poSrc->current.datfname.c_str());
     963           0 :         return CE_Failure;
     964             :     }
     965             : 
     966             :     // Write it then reissue the read
     967           1 :     err = poMRFDS->WriteTile(buf, infooffset, tinfo.size);
     968           1 :     CPLFree(buf);
     969           1 :     if (CE_None != err)
     970           0 :         return err;
     971             :     // Reissue read, it will work from the cloned data
     972           1 :     return IReadBlock(xblk, yblk, buffer);
     973             : }
     974             : 
     975             : /**
     976             :  *\brief read a block in the provided buffer
     977             :  *
     978             :  *  For separate band model, the DS buffer is not used, the read is direct in
     979             :  *the buffer For pixel interleaved model, the DS buffer holds the temp copy and
     980             :  *all the other bands are force read
     981             :  *
     982             :  */
     983             : 
     984        1877 : CPLErr MRFRasterBand::IReadBlock(int xblk, int yblk, void *buffer)
     985             : {
     986        1877 :     GInt32 cstride = img.pagesize.c;
     987             :     ILIdx tinfo;
     988        1877 :     ILSize req(xblk, yblk, 0, (nBand - 1) / cstride, m_l);
     989        3754 :     CPLDebug("MRF_IB",
     990             :              "IReadBlock %d,%d,0,%d, level %d, idxoffset " CPL_FRMT_GIB "\n",
     991        1877 :              xblk, yblk, nBand - 1, m_l, IdxOffset(req, img));
     992             : 
     993             :     // If this is a caching file and bypass is on, just do the fetch
     994        1877 :     if (poMRFDS->bypass_cache && !poMRFDS->source.empty())
     995           0 :         return FetchBlock(xblk, yblk, buffer);
     996             : 
     997        1877 :     tinfo.size = 0;  // Just in case it is missing
     998        1877 :     if (CE_None != poMRFDS->ReadTileIdx(tinfo, req, img))
     999             :     {
    1000           0 :         if (!poMRFDS->no_errors)
    1001             :         {
    1002           0 :             CPLError(CE_Failure, CPLE_AppDefined,
    1003             :                      "MRF: Unable to read index at offset " CPL_FRMT_GIB,
    1004           0 :                      IdxOffset(req, img));
    1005           0 :             return CE_Failure;
    1006             :         }
    1007           0 :         return FillBlock(buffer);
    1008             :     }
    1009             : 
    1010        1877 :     if (0 == tinfo.size)
    1011             :     {  // Could be missing or it could be caching
    1012             :         // Offset != 0 means no data, Update mode is for local MRFs only
    1013             :         // if caching index mode is RO don't try to fetch
    1014             :         // Also, caching MRFs can't be opened in update mode
    1015          12 :         if (0 != tinfo.offset || GA_Update == poMRFDS->eAccess ||
    1016          18 :             poMRFDS->source.empty() || IdxMode() == GF_Read)
    1017           2 :             return FillBlock(buffer);
    1018             : 
    1019             :         // caching MRF, need to fetch a block
    1020           4 :         return FetchBlock(xblk, yblk, buffer);
    1021             :     }
    1022             : 
    1023        1871 :     CPLDebug("MRF_IB", "Tinfo offset " CPL_FRMT_GIB ", size " CPL_FRMT_GIB "\n",
    1024             :              tinfo.offset, tinfo.size);
    1025             :     // If we have a tile, read it
    1026             : 
    1027             :     // Should use a permanent buffer, like the pbuffer mechanism
    1028             :     // Get a large buffer, in case we need to unzip
    1029             : 
    1030             :     // We add a padding of 3 bytes since in LERC1 decompression, we can
    1031             :     // dereference a unsigned int at the end of the buffer, that can be
    1032             :     // partially out of the buffer.
    1033             :     // Can be reproduced with :
    1034             :     // gdal_translate ../autotest/gcore/data/byte.tif out.mrf -of MRF -co
    1035             :     // COMPRESS=LERC -co OPTIONS=V1:YES -ot Float32 valgrind gdalinfo -checksum
    1036             :     // out.mrf Invalid read of size 4 at BitStuffer::read(unsigned char**,
    1037             :     // std::vector<unsigned int, std::allocator<unsigned int> >&) const
    1038             :     // (BitStuffer.cpp:153)
    1039             : 
    1040             :     // No stored tile should be larger than twice the raw size.
    1041        1871 :     if (tinfo.size <= 0 || tinfo.size > poMRFDS->pbsize * 2)
    1042             :     {
    1043           0 :         if (!poMRFDS->no_errors)
    1044             :         {
    1045           0 :             CPLError(CE_Failure, CPLE_OutOfMemory,
    1046             :                      "Stored tile is too large: " CPL_FRMT_GIB, tinfo.size);
    1047           0 :             return CE_Failure;
    1048             :         }
    1049           0 :         return FillBlock(buffer);
    1050             :     }
    1051             : 
    1052        1871 :     VSILFILE *dfp = DataFP();
    1053             : 
    1054             :     // No data file to read from
    1055        1871 :     if (dfp == nullptr)
    1056           0 :         return CE_Failure;
    1057             : 
    1058        1871 :     void *data = VSIMalloc(static_cast<size_t>(tinfo.size + PADDING_BYTES));
    1059        1871 :     if (data == nullptr)
    1060             :     {
    1061           0 :         CPLError(CE_Failure, CPLE_OutOfMemory,
    1062             :                  "Could not allocate memory for tile size: " CPL_FRMT_GIB,
    1063             :                  tinfo.size);
    1064           0 :         return CE_Failure;
    1065             :     }
    1066             : 
    1067             :     // This part is not thread safe, but it is what GDAL expects
    1068        1871 :     VSIFSeekL(dfp, tinfo.offset, SEEK_SET);
    1069        1871 :     if (1 != VSIFReadL(data, static_cast<size_t>(tinfo.size), 1, dfp))
    1070             :     {
    1071           0 :         CPLFree(data);
    1072           0 :         if (poMRFDS->no_errors)
    1073           0 :             return FillBlock(buffer);
    1074           0 :         CPLError(CE_Failure, CPLE_AppDefined, "Unable to read data page, %d@%x",
    1075           0 :                  static_cast<int>(tinfo.size), static_cast<int>(tinfo.offset));
    1076           0 :         return CE_Failure;
    1077             :     }
    1078             : 
    1079             :     /* initialize padding bytes */
    1080        1871 :     memset(((char *)data) + static_cast<size_t>(tinfo.size), 0, PADDING_BYTES);
    1081        1871 :     buf_mgr src = {(char *)data, static_cast<size_t>(tinfo.size)};
    1082             :     buf_mgr dst;
    1083             : 
    1084        1871 :     auto start_time = steady_clock::now();
    1085             : 
    1086             :     // We got the data, do we need to decompress it before decoding?
    1087        1871 :     if (dodeflate)
    1088             :     {
    1089           7 :         if (img.pageSizeBytes > INT_MAX - 1440)
    1090             :         {
    1091           0 :             CPLFree(data);
    1092           0 :             CPLError(CE_Failure, CPLE_AppDefined, "Page size is too big at %d",
    1093             :                      img.pageSizeBytes);
    1094           0 :             return CE_Failure;
    1095             :         }
    1096           7 :         dst.size =
    1097           7 :             img.pageSizeBytes +
    1098             :             1440;  // in case the packed page is a bit larger than the raw one
    1099           7 :         dst.buffer = static_cast<char *>(VSIMalloc(dst.size));
    1100           7 :         if (nullptr == dst.buffer)
    1101             :         {
    1102           0 :             CPLFree(data);
    1103           0 :             CPLError(CE_Failure, CPLE_OutOfMemory, "Cannot allocate %d bytes",
    1104           0 :                      static_cast<int>(dst.size));
    1105           0 :             return CE_Failure;
    1106             :         }
    1107             : 
    1108           7 :         if (ZUnPack(src, dst, deflate_flags))
    1109             :         {  // Got it unpacked, update the pointers
    1110           7 :             CPLFree(data);
    1111           7 :             data = dst.buffer;
    1112           7 :             tinfo.size = dst.size;
    1113             :         }
    1114             :         else
    1115             :         {  // assume the page was not gzipped, warn only
    1116           0 :             CPLFree(dst.buffer);
    1117           0 :             if (!poMRFDS->no_errors)
    1118           0 :                 CPLError(CE_Warning, CPLE_AppDefined, "Can't inflate page!");
    1119             :         }
    1120             :     }
    1121             : 
    1122             : #if defined(ZSTD_SUPPORT)
    1123             :     // undo ZSTD
    1124        1864 :     else if (dozstd)
    1125             :     {
    1126           9 :         auto ctx = poMRFDS->getzsd();
    1127           9 :         if (!ctx)
    1128             :         {
    1129           0 :             CPLFree(data);
    1130           0 :             CPLError(CE_Failure, CPLE_AppDefined, "Can't acquire ZSTD context");
    1131           0 :             return CE_Failure;
    1132             :         }
    1133           9 :         if (img.pageSizeBytes > INT_MAX - 1440)
    1134             :         {
    1135           0 :             CPLFree(data);
    1136           0 :             CPLError(CE_Failure, CPLE_AppDefined, "Page is too large at %d",
    1137             :                      img.pageSizeBytes);
    1138           0 :             return CE_Failure;
    1139             :         }
    1140           9 :         dst.size =
    1141           9 :             img.pageSizeBytes +
    1142             :             1440;  // Allow for a slight increase from previous compressions
    1143           9 :         dst.buffer = static_cast<char *>(VSIMalloc(dst.size));
    1144           9 :         if (nullptr == dst.buffer)
    1145             :         {
    1146           0 :             CPLFree(data);
    1147           0 :             CPLError(CE_Failure, CPLE_OutOfMemory, "Cannot allocate %d bytes",
    1148           0 :                      static_cast<int>(dst.size));
    1149           0 :             return CE_Failure;
    1150             :         }
    1151             : 
    1152          18 :         auto raw_size = ZSTD_decompressDCtx(ctx, dst.buffer, dst.size,
    1153           9 :                                             src.buffer, src.size);
    1154           9 :         if (ZSTD_isError(raw_size))
    1155             :         {  // assume page was not packed, warn only
    1156           0 :             CPLFree(dst.buffer);
    1157           0 :             if (!poMRFDS->no_errors)
    1158           0 :                 CPLError(CE_Warning, CPLE_AppDefined,
    1159             :                          "Can't unpack ZSTD page!");
    1160             :         }
    1161             :         else
    1162             :         {
    1163           9 :             CPLFree(data);  // The compressed data
    1164           9 :             data = dst.buffer;
    1165           9 :             tinfo.size = raw_size;
    1166             :             // Might need to undo the rank sort
    1167           9 :             size_t ranks = 0;
    1168           9 :             if (img.comp == IL_NONE || img.comp == IL_ZSTD)
    1169           9 :                 ranks = static_cast<size_t>(GDALGetDataTypeSizeBytes(img.dt)) *
    1170           9 :                         img.pagesize.c;
    1171           9 :             if (ranks)
    1172             :             {
    1173           9 :                 src.buffer = static_cast<char *>(data);
    1174           9 :                 src.size = static_cast<size_t>(tinfo.size);
    1175           9 :                 derank(src, ranks);
    1176             :             }
    1177             :         }
    1178             :     }
    1179             : #endif
    1180             : 
    1181        1871 :     src.buffer = static_cast<char *>(data);
    1182        1871 :     src.size = static_cast<size_t>(tinfo.size);
    1183             : 
    1184             :     // After unpacking, the size has to be pageSizeBytes
    1185             :     // If pages are interleaved, use the dataset page buffer instead
    1186        1871 :     dst.buffer = reinterpret_cast<char *>(
    1187           8 :         (1 == cstride) ? buffer : poMRFDS->GetPBuffer());
    1188        1871 :     dst.size = img.pageSizeBytes;
    1189             : 
    1190        1871 :     if (poMRFDS->no_errors)
    1191           0 :         CPLPushErrorHandler(CPLQuietErrorHandler);
    1192        1871 :     CPLErr ret = Decompress(dst, src);
    1193             : 
    1194        1871 :     poMRFDS->read_timer +=
    1195        1871 :         duration_cast<nanoseconds>(steady_clock::now() - start_time);
    1196             : 
    1197        1871 :     dst.size =
    1198        1871 :         img.pageSizeBytes;  // In case the decompress failed, force it back
    1199             : 
    1200             :     // Swap whatever we decompressed if we need to
    1201        1871 :     if (is_Endianness_Dependent(img.dt, img.comp) && (img.nbo != NET_ORDER))
    1202           0 :         swab_buff(dst, img);
    1203             : 
    1204        1871 :     CPLFree(data);
    1205        1871 :     if (poMRFDS->no_errors)
    1206             :     {
    1207           0 :         CPLPopErrorHandler();
    1208           0 :         if (ret != CE_None)  // Set each page buffer to the correct no data
    1209             :                              // value, then proceed
    1210           0 :             return (1 == cstride) ? FillBlock(buffer)
    1211           0 :                                   : FillBlock(xblk, yblk, buffer);
    1212             :     }
    1213             : 
    1214             :     // If pages are separate or we had errors, we're done
    1215        1871 :     if (1 == cstride || CE_None != ret)
    1216        1863 :         return ret;
    1217             : 
    1218             :     // De-interleave page from dataset buffer and return
    1219           8 :     return ReadInterleavedBlock(xblk, yblk, buffer);
    1220             : }
    1221             : 
    1222             : /**
    1223             :  *\brief Write a block from the provided buffer
    1224             :  *
    1225             :  * Same trick as read, use a temporary tile buffer for pixel interleave
    1226             :  * For band separate, use a
    1227             :  * Write the block once it has all the bands, report
    1228             :  * if a new block is started before the old one was completed
    1229             :  *
    1230             :  */
    1231             : 
    1232        5076 : CPLErr MRFRasterBand::IWriteBlock(int xblk, int yblk, void *buffer)
    1233             : {
    1234        5076 :     GInt32 cstride = img.pagesize.c;
    1235        5076 :     ILSize req(xblk, yblk, 0, (nBand - 1) / cstride, m_l);
    1236        5076 :     GUIntBig infooffset = IdxOffset(req, img);
    1237             : 
    1238        5076 :     CPLDebug("MRF_IB", "IWriteBlock %d,%d,0,%d, level %d, stride %d\n", xblk,
    1239             :              yblk, nBand, m_l, cstride);
    1240             : 
    1241             :     // Finish the Create call
    1242        5076 :     if (!poMRFDS->bCrystalized && !poMRFDS->Crystalize())
    1243             :     {
    1244           0 :         CPLError(CE_Failure, CPLE_AppDefined, "MRF: Error creating files");
    1245           0 :         return CE_Failure;
    1246             :     }
    1247             : 
    1248        5076 :     if (1 == cstride)
    1249             :     {  // Separate bands, we can write it as is
    1250             :         // Empty page skip
    1251             :         int success;
    1252        5065 :         double val = GetNoDataValue(&success);
    1253        5065 :         if (!success)
    1254        4993 :             val = 0.0;
    1255        5065 :         if (isAllVal(eDataType, buffer, img.pageSizeBytes, val))
    1256           1 :             return poMRFDS->WriteTile(nullptr, infooffset, 0);
    1257             : 
    1258             :         // Use the pbuffer to hold the compressed page before writing it
    1259        5064 :         poMRFDS->tile = ILSize();  // Mark it corrupt
    1260             : 
    1261             :         buf_mgr src;
    1262        5064 :         src.buffer = (char *)buffer;
    1263        5064 :         src.size = static_cast<size_t>(img.pageSizeBytes);
    1264        5064 :         buf_mgr dst = {(char *)poMRFDS->GetPBuffer(),
    1265        5064 :                        poMRFDS->GetPBufferSize()};
    1266             : 
    1267             :         // Swab the source before encoding if we need to
    1268        5064 :         if (is_Endianness_Dependent(img.dt, img.comp) && (img.nbo != NET_ORDER))
    1269           0 :             swab_buff(src, img);
    1270             : 
    1271        5064 :         auto start_time = steady_clock::now();
    1272             : 
    1273             :         // Compress functions need to return the compressed size in
    1274             :         // the bytes in buffer field
    1275        5064 :         Compress(dst, src);
    1276        5064 :         void *usebuff = dst.buffer;
    1277        5064 :         if (dodeflate)
    1278             :         {
    1279             :             usebuff =
    1280           7 :                 DeflateBlock(dst, poMRFDS->pbsize - dst.size, deflate_flags);
    1281           7 :             if (!usebuff)
    1282             :             {
    1283           0 :                 CPLError(CE_Failure, CPLE_AppDefined, "MRF: Deflate error");
    1284           0 :                 return CE_Failure;
    1285             :             }
    1286             :         }
    1287             : 
    1288             : #if defined(ZSTD_SUPPORT)
    1289        5057 :         else if (dozstd)
    1290             :         {
    1291          10 :             size_t ranks = 0;  // Assume no need for byte rank sort
    1292          10 :             if (img.comp == IL_NONE || img.comp == IL_ZSTD)
    1293          10 :                 ranks = static_cast<size_t>(GDALGetDataTypeSizeBytes(img.dt));
    1294          10 :             usebuff = ZstdCompBlock(dst, poMRFDS->pbsize - dst.size, zstd_level,
    1295          10 :                                     poMRFDS->getzsc(), ranks);
    1296          10 :             if (!usebuff)
    1297             :             {
    1298           0 :                 CPLError(CE_Failure, CPLE_AppDefined,
    1299             :                          "MRF: Zstd Compression error");
    1300           0 :                 return CE_Failure;
    1301             :             }
    1302             :         }
    1303             : #endif
    1304        5064 :         poMRFDS->write_timer +=
    1305        5064 :             duration_cast<nanoseconds>(steady_clock::now() - start_time);
    1306        5064 :         return poMRFDS->WriteTile(usebuff, infooffset, dst.size);
    1307             :     }
    1308             : 
    1309             :     // Multiple bands per page, use a temporary to assemble the page
    1310             :     // Temporary is large because we use it to hold both the uncompressed and
    1311             :     // the compressed
    1312          11 :     poMRFDS->tile = req;
    1313          11 :     poMRFDS->bdirty = 0;
    1314             : 
    1315             :     // Keep track of what bands are empty
    1316          11 :     GUIntBig empties = 0;
    1317             : 
    1318          11 :     void *tbuffer = VSIMalloc(img.pageSizeBytes + poMRFDS->pbsize);
    1319             : 
    1320          11 :     if (!tbuffer)
    1321             :     {
    1322           0 :         CPLError(CE_Failure, CPLE_AppDefined,
    1323             :                  "MRF: Can't allocate write buffer");
    1324           0 :         return CE_Failure;
    1325             :     }
    1326             : 
    1327             :     // Get the other bands from the block cache
    1328          44 :     for (int iBand = 0; iBand < poMRFDS->nBands; iBand++)
    1329             :     {
    1330          33 :         char *pabyThisImage = nullptr;
    1331          33 :         GDALRasterBlock *poBlock = nullptr;
    1332             : 
    1333          33 :         if (iBand == nBand - 1)
    1334             :         {
    1335          11 :             pabyThisImage = reinterpret_cast<char *>(buffer);
    1336          11 :             poMRFDS->bdirty |= bandbit();
    1337             :         }
    1338             :         else
    1339             :         {
    1340          22 :             GDALRasterBand *band = poMRFDS->GetRasterBand(iBand + 1);
    1341             :             // Pick the right overview
    1342          22 :             if (m_l)
    1343           0 :                 band = band->GetOverview(m_l - 1);
    1344             :             poBlock = (reinterpret_cast<MRFRasterBand *>(band))
    1345          22 :                           ->TryGetLockedBlockRef(xblk, yblk);
    1346          22 :             if (nullptr == poBlock)
    1347           0 :                 continue;
    1348             :             // This is where the image data is for this band
    1349             : 
    1350          22 :             pabyThisImage = reinterpret_cast<char *>(poBlock->GetDataRef());
    1351          22 :             poMRFDS->bdirty |= bandbit(iBand);
    1352             :         }
    1353             : 
    1354             :         // Keep track of empty bands, but encode them anyhow, in case some are
    1355             :         // not empty
    1356             :         int success;
    1357          33 :         double val = GetNoDataValue(&success);
    1358          33 :         if (!success)
    1359          33 :             val = 0.0;
    1360          33 :         if (isAllVal(eDataType, pabyThisImage, blockSizeBytes(), val))
    1361           6 :             empties |= bandbit(iBand);
    1362             : 
    1363             :             // Copy the data into the dataset buffer here
    1364             :             // Just the right mix of templates and macros make this real tidy
    1365             : #define CpySO(T)                                                               \
    1366             :     cpy_stride_out<T>((reinterpret_cast<T *>(tbuffer)) + iBand, pabyThisImage, \
    1367             :                       blockSizeBytes() / sizeof(T), cstride)
    1368             : 
    1369             :         // Build the page in tbuffer
    1370          33 :         switch (GDALGetDataTypeSize(eDataType) / 8)
    1371             :         {
    1372          33 :             case 1:
    1373          33 :                 CpySO(GByte);
    1374          33 :                 break;
    1375           0 :             case 2:
    1376           0 :                 CpySO(GInt16);
    1377           0 :                 break;
    1378           0 :             case 4:
    1379           0 :                 CpySO(GInt32);
    1380           0 :                 break;
    1381           0 :             case 8:
    1382           0 :                 CpySO(GIntBig);
    1383           0 :                 break;
    1384           0 :             default:
    1385             :             {
    1386           0 :                 CPLError(CE_Failure, CPLE_AppDefined,
    1387             :                          "MRF: Write datatype of %d bytes "
    1388             :                          "not implemented",
    1389           0 :                          GDALGetDataTypeSize(eDataType) / 8);
    1390           0 :                 if (poBlock != nullptr)
    1391             :                 {
    1392           0 :                     poBlock->MarkClean();
    1393           0 :                     poBlock->DropLock();
    1394             :                 }
    1395           0 :                 CPLFree(tbuffer);
    1396           0 :                 return CE_Failure;
    1397             :             }
    1398             :         }
    1399             : 
    1400          33 :         if (poBlock != nullptr)
    1401             :         {
    1402          22 :             poBlock->MarkClean();
    1403          22 :             poBlock->DropLock();
    1404             :         }
    1405             :     }
    1406             : 
    1407             :     // Should keep track of the individual band buffers and only mix them if
    1408             :     // this is not an empty page ( move the Copy with Stride Out from above
    1409             :     // below this test This way works fine, but it does work extra for empty
    1410             :     // pages
    1411             : 
    1412          11 :     if (GIntBig(empties) == AllBandMask())
    1413             :     {
    1414           0 :         CPLFree(tbuffer);
    1415           0 :         return poMRFDS->WriteTile(nullptr, infooffset, 0);
    1416             :     }
    1417             : 
    1418          11 :     if (poMRFDS->bdirty != AllBandMask())
    1419           0 :         CPLError(CE_Warning, CPLE_AppDefined,
    1420             :                  "MRF: IWrite, band dirty mask is " CPL_FRMT_GIB
    1421             :                  " instead of " CPL_FRMT_GIB,
    1422           0 :                  poMRFDS->bdirty, AllBandMask());
    1423             : 
    1424             :     buf_mgr src;
    1425          11 :     src.buffer = (char *)tbuffer;
    1426          11 :     src.size = static_cast<size_t>(img.pageSizeBytes);
    1427             : 
    1428             :     // Use the space after pagesizebytes for compressed output, it is of pbsize
    1429          11 :     char *outbuff = (char *)tbuffer + img.pageSizeBytes;
    1430             : 
    1431          11 :     buf_mgr dst = {outbuff, poMRFDS->pbsize};
    1432             :     CPLErr ret;
    1433             : 
    1434          11 :     auto start_time = steady_clock::now();
    1435             : 
    1436          11 :     ret = Compress(dst, src);
    1437          11 :     if (ret != CE_None)
    1438             :     {
    1439             :         // Compress failed, write it as an empty tile
    1440           0 :         CPLFree(tbuffer);
    1441           0 :         poMRFDS->WriteTile(nullptr, infooffset, 0);
    1442           0 :         return CE_None;  // Should report the error, but it triggers partial
    1443             :                          // band attempts
    1444             :     }
    1445             : 
    1446             :     // Where the output is, in case we deflate
    1447          11 :     void *usebuff = outbuff;
    1448          11 :     if (dodeflate)
    1449             :     {
    1450             :         // Move the packed part at the start of tbuffer, to make more space
    1451             :         // available
    1452           0 :         memcpy(tbuffer, outbuff, dst.size);
    1453           0 :         dst.buffer = static_cast<char *>(tbuffer);
    1454           0 :         usebuff = DeflateBlock(dst,
    1455           0 :                                static_cast<size_t>(img.pageSizeBytes) +
    1456           0 :                                    poMRFDS->pbsize - dst.size,
    1457             :                                deflate_flags);
    1458           0 :         if (!usebuff)
    1459           0 :             CPLError(CE_Failure, CPLE_AppDefined, "MRF: Deflate error");
    1460             :     }
    1461             : 
    1462             : #if defined(ZSTD_SUPPORT)
    1463          11 :     else if (dozstd)
    1464             :     {
    1465           1 :         memcpy(tbuffer, outbuff, dst.size);
    1466           1 :         dst.buffer = static_cast<char *>(tbuffer);
    1467           1 :         size_t ranks = 0;  // Assume no need for byte rank sort
    1468           1 :         if (img.comp == IL_NONE || img.comp == IL_ZSTD)
    1469           1 :             ranks =
    1470           1 :                 static_cast<size_t>(GDALGetDataTypeSizeBytes(img.dt)) * cstride;
    1471           3 :         usebuff = ZstdCompBlock(dst,
    1472           1 :                                 static_cast<size_t>(img.pageSizeBytes) +
    1473           1 :                                     poMRFDS->pbsize - dst.size,
    1474           1 :                                 zstd_level, poMRFDS->getzsc(), ranks);
    1475           1 :         if (!usebuff)
    1476           0 :             CPLError(CE_Failure, CPLE_AppDefined,
    1477             :                      "MRF: ZStd compression error");
    1478             :     }
    1479             : #endif
    1480             : 
    1481          11 :     poMRFDS->write_timer +=
    1482          11 :         duration_cast<nanoseconds>(steady_clock::now() - start_time);
    1483             : 
    1484          11 :     if (!usebuff)
    1485             :     {  // Error was signaled
    1486           0 :         CPLFree(tbuffer);
    1487           0 :         poMRFDS->WriteTile(nullptr, infooffset, 0);
    1488           0 :         poMRFDS->bdirty = 0;
    1489           0 :         return CE_Failure;
    1490             :     }
    1491             : 
    1492          11 :     ret = poMRFDS->WriteTile(usebuff, infooffset, dst.size);
    1493          11 :     CPLFree(tbuffer);
    1494             : 
    1495          11 :     poMRFDS->bdirty = 0;
    1496          11 :     return ret;
    1497             : }
    1498             : 
    1499             : //
    1500             : // Tests if a given block exists without reading it
    1501             : // returns false only when it is definitely not existing
    1502             : //
    1503         160 : bool MRFRasterBand::TestBlock(int xblk, int yblk)
    1504             : {
    1505             :     // When bypassing the cache, assume all blocks are valid
    1506         160 :     if (poMRFDS->bypass_cache && !poMRFDS->source.empty())
    1507           0 :         return true;
    1508             : 
    1509             :     // Blocks outside of image have no data by default
    1510         160 :     if (xblk < 0 || yblk < 0 || xblk >= img.pagecount.x ||
    1511         146 :         yblk >= img.pagecount.y)
    1512          25 :         return false;
    1513             : 
    1514             :     ILIdx tinfo;
    1515         135 :     GInt32 cstride = img.pagesize.c;
    1516         135 :     ILSize req(xblk, yblk, 0, (nBand - 1) / cstride, m_l);
    1517             : 
    1518         135 :     if (CE_None != poMRFDS->ReadTileIdx(tinfo, req, img))
    1519             :         // Got an error reading the tile index
    1520           0 :         return !poMRFDS->no_errors;
    1521             : 
    1522             :     // Got an index, if the size is readable, the block does exist
    1523         135 :     if (0 < tinfo.size && tinfo.size < poMRFDS->pbsize * 2)
    1524         135 :         return true;
    1525             : 
    1526             :     // We are caching, but the tile has not been checked, so it could exist
    1527           0 :     return (!poMRFDS->source.empty() && 0 == tinfo.offset);
    1528             : }
    1529             : 
    1530         183 : int MRFRasterBand::GetOverviewCount()
    1531             : {
    1532             :     // First try internal overviews
    1533         183 :     int nInternalOverviewCount = static_cast<int>(overviews.size());
    1534         183 :     if (nInternalOverviewCount > 0)
    1535         137 :         return nInternalOverviewCount;
    1536          46 :     return GDALPamRasterBand::GetOverviewCount();
    1537             : }
    1538             : 
    1539         158 : GDALRasterBand *MRFRasterBand::GetOverview(int n)
    1540             : {
    1541             :     // First try internal overviews
    1542         158 :     if (n >= 0 && n < (int)overviews.size())
    1543          85 :         return overviews[n];
    1544          73 :     return GDALPamRasterBand::GetOverview(n);
    1545             : }
    1546             : 
    1547             : NAMESPACE_MRF_END

Generated by: LCOV version 1.14