1 /*
2 RawSpeed - RAW file decoder.
3
4 Copyright (C) 2018 Stefan Löffler
5 Copyright (C) 2018-2019 Roman Lebedev
6
7 This library is free software; you can redistribute it and/or
8 modify it under the terms of the GNU Lesser General Public
9 License as published by the Free Software Foundation; either
10 version 2 of the License, or (at your option) any later version.
11
12 This library is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
15 Lesser General Public License for more details.
16
17 You should have received a copy of the GNU Lesser General Public
18 License along with this library; if not, write to the Free Software
19 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
20 */
21
22 /*
23 This is a decompressor for VC-5 raw compression algo, as used in GoPro raws.
24 This implementation is similar to that one of the official reference
25 implementation of the https://github.com/gopro/gpr project, and is producing
26 bitwise-identical output as compared with the Adobe DNG Converter
27 implementation.
28 */
29
30 #include "rawspeedconfig.h" // for HAVE_OPENMP
31 #include "decompressors/VC5Decompressor.h"
32 #include "common/Array2DRef.h" // for Array2DRef
33 #include "common/Common.h" // for clampBits, roundUpDivision
34 #include "common/Point.h" // for iPoint2D
35 #include "common/RawspeedException.h" // for RawspeedException
36 #include "common/SimpleLUT.h" // for SimpleLUT, SimpleLUT<>::va...
37 #include "decoders/RawDecoderException.h" // for ThrowRDE
38 #include "io/Endianness.h" // for Endianness, Endianness::big
39 #include <cassert> // for assert
40 #include <cmath> // for pow
41 #include <initializer_list> // for initializer_list
42 #include <limits> // for numeric_limits
43 #include <optional> // for optional
44 #include <string> // for string
45 #include <utility> // for move
46
47 namespace {
48
49 // Definitions needed by table17.inc
50 // Taken from
51 // https://github.com/gopro/gpr/blob/a513701afce7b03173213a2f67dfd9dd28fa1868/source/lib/vc5_decoder/vlc.h
52 struct RLV {
53 uint_fast8_t size; //!< Size of code word in bits
54 uint32_t bits; //!< Code word bits right justified
55 uint16_t count; //!< Run length
56 uint16_t value; //!< Run value (unsigned)
57 };
58 #define RLVTABLE(n) \
59 struct { \
60 const uint32_t length; \
61 const RLV entries[n]; \
62 } constexpr
63 #include "gopro/vc5/table17.inc"
64
decompand(int16_t val)65 constexpr int16_t decompand(int16_t val) {
66 double c = val;
67 // Invert companding curve
68 c += (c * c * c * 768) / (255. * 255. * 255.);
69 if (c > std::numeric_limits<int16_t>::max())
70 return std::numeric_limits<int16_t>::max();
71 if (c < std::numeric_limits<int16_t>::min())
72 return std::numeric_limits<int16_t>::min();
73 return c;
74 }
75
76 #ifndef NDEBUG
__anon8f7636770202() 77 [[maybe_unused]] const int ignore = []() {
78 for (const RLV& entry : table17.entries) {
79 assert(((-decompand(entry.value)) == decompand(-int16_t(entry.value))) &&
80 "negation of decompanded value is the same as decompanding of "
81 "negated value");
82 }
83 return 0;
84 }();
85 #endif
86
__anon8f7636770302() 87 const std::array<RLV, table17.length> decompandedTable17 = []() {
88 std::array<RLV, table17.length> d;
89 for (auto i = 0U; i < table17.length; i++) {
90 d[i] = table17.entries[i];
91 d[i].value = decompand(table17.entries[i].value);
92 }
93 return d;
94 }();
95
readValue(bool & storage)96 inline bool readValue(bool& storage) {
97 bool value;
98
99 #ifdef HAVE_OPENMP
100 #pragma omp atomic read
101 #endif
102 value = storage;
103
104 return value;
105 }
106
107 } // namespace
108
109 #define PRECISION_MIN 8
110 #define PRECISION_MAX 16
111
112 #define MARKER_BAND_END 1
113
114 namespace rawspeed {
115
setBandValid(const int band)116 void VC5Decompressor::Wavelet::setBandValid(const int band) {
117 mDecodedBandMask |= (1 << band);
118 }
119
isBandValid(const int band) const120 bool VC5Decompressor::Wavelet::isBandValid(const int band) const {
121 return mDecodedBandMask & (1 << band);
122 }
123
allBandsValid() const124 bool VC5Decompressor::Wavelet::allBandsValid() const {
125 return mDecodedBandMask == static_cast<uint32_t>((1 << maxBands) - 1);
126 }
127
128 namespace {
129 template <typename LowGetter>
convolute(int row,int col,std::array<int,4> muls,const Array2DRef<const int16_t> high,LowGetter lowGetter,int DescaleShift=0)130 static inline auto convolute(int row, int col, std::array<int, 4> muls,
131 const Array2DRef<const int16_t> high,
132 LowGetter lowGetter, int DescaleShift = 0) {
133 auto highCombined = muls[0] * high(row, col);
134 auto lowsCombined = [muls, lowGetter]() {
135 int lows = 0;
136 for (int i = 0; i < 3; i++)
137 lows += muls[1 + i] * lowGetter(i);
138 return lows;
139 }();
140 // Round up 'lows' up
141 lowsCombined += 4;
142 // And finally 'average' them.
143 auto lowsRounded = lowsCombined >> 3;
144 auto total = highCombined + lowsRounded;
145 // Descale it.
146 total <<= DescaleShift;
147 // And average it.
148 total >>= 1;
149 return total;
150 }
151
152 struct ConvolutionParams {
153 struct First {
154 static constexpr std::array<int, 4> mul_even = {+1, +11, -4, +1};
155 static constexpr std::array<int, 4> mul_odd = {-1, +5, +4, -1};
156 static constexpr int coord_shift = 0;
157 };
158
159 struct Middle {
160 static constexpr std::array<int, 4> mul_even = {+1, +1, +8, -1};
161 static constexpr std::array<int, 4> mul_odd = {-1, -1, +8, +1};
162 static constexpr int coord_shift = -1;
163 };
164
165 struct Last {
166 static constexpr std::array<int, 4> mul_even = {+1, -1, +4, +5};
167 static constexpr std::array<int, 4> mul_odd = {-1, +1, -4, +11};
168 static constexpr int coord_shift = -2;
169 };
170 };
171
172 } // namespace
173
reconstructPass(const Array2DRef<const int16_t> high,const Array2DRef<const int16_t> low)174 VC5Decompressor::BandData VC5Decompressor::Wavelet::reconstructPass(
175 const Array2DRef<const int16_t> high,
176 const Array2DRef<const int16_t> low) noexcept {
177 BandData combined;
178 auto& dst = combined.description;
179 dst = Array2DRef<int16_t>::create(combined.storage, high.width,
180 2 * high.height);
181
182 auto process = [low, high, dst](auto segment, int row, int col) {
183 using SegmentTy = decltype(segment);
184 auto lowGetter = [&row, &col, low](int delta) {
185 return low(row + SegmentTy::coord_shift + delta, col);
186 };
187 auto convolution = [&row, &col, high, lowGetter](std::array<int, 4> muls) {
188 return convolute(row, col, muls, high, lowGetter, /*DescaleShift*/ 0);
189 };
190
191 int even = convolution(SegmentTy::mul_even);
192 int odd = convolution(SegmentTy::mul_odd);
193
194 dst(2 * row, col) = static_cast<int16_t>(even);
195 dst(2 * row + 1, col) = static_cast<int16_t>(odd);
196 };
197
198 #pragma GCC diagnostic push
199 // See https://bugs.llvm.org/show_bug.cgi?id=51666
200 #pragma GCC diagnostic ignored "-Wsign-compare"
201 #ifdef HAVE_OPENMP
202 #pragma omp taskloop default(none) firstprivate(dst, process) num_tasks( \
203 roundUpDivision(rawspeed_get_number_of_processor_cores(), numChannels))
204 #endif
205 for (int row = 0; row < dst.height / 2; ++row) {
206 #pragma GCC diagnostic pop
207 if (row == 0) {
208 // 1st row
209 for (int col = 0; col < dst.width; ++col)
210 process(ConvolutionParams::First(), row, col);
211 } else if (row + 1 < dst.height / 2) {
212 // middle rows
213 for (int col = 0; col < dst.width; ++col)
214 process(ConvolutionParams::Middle(), row, col);
215 } else {
216 // last row
217 for (int col = 0; col < dst.width; ++col)
218 process(ConvolutionParams::Last(), row, col);
219 }
220 }
221
222 return combined;
223 }
224
combineLowHighPass(const Array2DRef<const int16_t> low,const Array2DRef<const int16_t> high,int descaleShift,bool clampUint=false,bool finalWavelet=false)225 VC5Decompressor::BandData VC5Decompressor::Wavelet::combineLowHighPass(
226 const Array2DRef<const int16_t> low, const Array2DRef<const int16_t> high,
227 int descaleShift, bool clampUint = false,
228 bool finalWavelet = false) noexcept {
229 BandData combined;
230 auto& dst = combined.description;
231 dst = Array2DRef<int16_t>::create(combined.storage, 2 * high.width,
232 high.height);
233
234 auto process = [low, high, descaleShift, clampUint, dst](auto segment,
235 int row, int col) {
236 using SegmentTy = decltype(segment);
237 auto lowGetter = [&row, &col, low](int delta) {
238 return low(row, col + SegmentTy::coord_shift + delta);
239 };
240 auto convolution = [&row, &col, high, lowGetter,
241 descaleShift](std::array<int, 4> muls) {
242 return convolute(row, col, muls, high, lowGetter, descaleShift);
243 };
244
245 int even = convolution(SegmentTy::mul_even);
246 int odd = convolution(SegmentTy::mul_odd);
247
248 if (clampUint) {
249 even = clampBits(even, 14);
250 odd = clampBits(odd, 14);
251 }
252 dst(row, 2 * col) = static_cast<int16_t>(even);
253 dst(row, 2 * col + 1) = static_cast<int16_t>(odd);
254 };
255
256 // Horizontal reconstruction
257 #pragma GCC diagnostic push
258 // See https://bugs.llvm.org/show_bug.cgi?id=51666
259 #pragma GCC diagnostic ignored "-Wsign-compare"
260 #ifdef HAVE_OPENMP
261 #pragma omp taskloop if (finalWavelet) default(none) firstprivate(dst, \
262 process) \
263 num_tasks(roundUpDivision(rawspeed_get_number_of_processor_cores(), 2)) \
264 mergeable
265 #endif
266 for (int row = 0; row < dst.height; ++row) {
267 #pragma GCC diagnostic pop
268 // First col
269 int col = 0;
270 process(ConvolutionParams::First(), row, col);
271 // middle cols
272 for (col = 1; col + 1 < dst.width / 2; ++col) {
273 process(ConvolutionParams::Middle(), row, col);
274 }
275 // last col
276 process(ConvolutionParams::Last(), row, col);
277 }
278
279 return combined;
280 }
281
282 void VC5Decompressor::Wavelet::ReconstructableBand::
createLowpassReconstructionTask(bool & exceptionThrown)283 createLowpassReconstructionTask(bool& exceptionThrown) noexcept {
284 auto& highlow = wavelet.bands[2]->data;
285 auto& lowlow = wavelet.bands[0]->data;
286 auto& lowpass = intermediates.lowpass;
287
288 #ifdef HAVE_OPENMP
289 #pragma omp task default(none) \
290 shared(exceptionThrown, highlow, lowlow, lowpass) \
291 depend(in \
292 : highlow, lowlow) depend(out \
293 : lowpass)
294 #endif
295 {
296 // Proceed only if decoding did not fail.
297 if (!readValue(exceptionThrown)) {
298 assert(highlow.has_value() && lowlow.has_value() &&
299 "Failed to produce precursor bands?");
300 // Reconstruct the "immediates", the actual low pass ...
301 assert(!lowpass.has_value() && "Combined this precursor band already?");
302 lowpass.emplace(
303 Wavelet::reconstructPass(highlow->description, lowlow->description));
304 }
305 }
306 }
307
308 void VC5Decompressor::Wavelet::ReconstructableBand::
createHighpassReconstructionTask(bool & exceptionThrown)309 createHighpassReconstructionTask(bool& exceptionThrown) noexcept {
310 auto& highhigh = wavelet.bands[3]->data;
311 auto& lowhigh = wavelet.bands[1]->data;
312 auto& highpass = intermediates.highpass;
313
314 #ifdef HAVE_OPENMP
315 #pragma omp task default(none) \
316 shared(exceptionThrown, highhigh, lowhigh, highpass) \
317 depend(in \
318 : highhigh, lowhigh) depend(out \
319 : highpass)
320 #endif
321 {
322 // Proceed only if decoding did not fail.
323 if (!readValue(exceptionThrown)) {
324 assert(highhigh.has_value() && lowhigh.has_value() &&
325 "Failed to produce precursor bands?");
326 assert(!highpass.has_value() && "Combined this precursor band already?");
327 highpass.emplace(Wavelet::reconstructPass(highhigh->description,
328 lowhigh->description));
329 }
330 }
331 }
332
333 void VC5Decompressor::Wavelet::ReconstructableBand::
createLowHighPassCombiningTask(bool & exceptionThrown)334 createLowHighPassCombiningTask(bool& exceptionThrown) noexcept {
335 auto& lowpass = intermediates.lowpass;
336 auto& highpass = intermediates.highpass;
337 auto& reconstructedLowpass = data;
338
339 #ifdef HAVE_OPENMP
340 #pragma omp task default(none) depend(in : lowpass, highpass)
341 #endif
342 wavelet.bands.clear();
343
344 #ifdef HAVE_OPENMP
345 #pragma omp task default(none) \
346 shared(exceptionThrown, lowpass, highpass, reconstructedLowpass) \
347 depend(in \
348 : lowpass, highpass) depend(out \
349 : reconstructedLowpass)
350 #endif
351 {
352 // Proceed only if decoding did not fail.
353 if (!readValue(exceptionThrown)) {
354 assert(lowpass.has_value() && highpass.has_value() &&
355 "Failed to combine precursor bands?");
356 assert(!data.has_value() && "Reconstructed this band already?");
357
358 int16_t descaleShift = (wavelet.prescale == 2 ? 2 : 0);
359
360 // And finally, combine the low pass, and high pass.
361 reconstructedLowpass.emplace(Wavelet::combineLowHighPass(
362 lowpass->description, highpass->description, descaleShift, clampUint,
363 finalWavelet));
364 }
365 }
366 }
367
createDecodingTasks(ErrorLog & errLog,bool & exceptionThrow)368 void VC5Decompressor::Wavelet::ReconstructableBand::createDecodingTasks(
369 ErrorLog& errLog, bool& exceptionThrow) noexcept {
370 assert(wavelet.allBandsValid());
371 createLowpassReconstructionTask(exceptionThrow);
372 createHighpassReconstructionTask(exceptionThrow);
373 createLowHighPassCombiningTask(exceptionThrow);
374 }
375
VC5Decompressor(ByteStream bs,const RawImage & img)376 VC5Decompressor::VC5Decompressor(ByteStream bs, const RawImage& img)
377 : mRaw(img), mBs(std::move(bs)) {
378 if (!mRaw->dim.hasPositiveArea())
379 ThrowRDE("Bad image dimensions.");
380
381 if (mRaw->dim.x % mVC5.patternWidth != 0)
382 ThrowRDE("Width %u is not a multiple of %u", mRaw->dim.x,
383 mVC5.patternWidth);
384
385 if (mRaw->dim.y % mVC5.patternHeight != 0)
386 ThrowRDE("Height %u is not a multiple of %u", mRaw->dim.y,
387 mVC5.patternHeight);
388
389 // Initialize wavelet sizes.
390 for (Channel& channel : channels) {
391 uint16_t waveletWidth = mRaw->dim.x;
392 uint16_t waveletHeight = mRaw->dim.y;
393 for (Wavelet& wavelet : channel.wavelets) {
394 // Pad dimensions as necessary and divide them by two for the next wavelet
395 for (auto* dimension : {&waveletWidth, &waveletHeight})
396 *dimension = roundUpDivision(*dimension, 2);
397 wavelet.width = waveletWidth;
398 wavelet.height = waveletHeight;
399
400 wavelet.bands.resize(
401 &wavelet == channel.wavelets.begin() ? 1 : Wavelet::maxBands);
402 }
403 }
404
405 if (img->whitePoint <= 0 || img->whitePoint > int(((1U << 16U) - 1U)))
406 ThrowRDE("Bad white level %i", img->whitePoint);
407
408 outputBits = 0;
409 for (int wp = img->whitePoint; wp != 0; wp >>= 1)
410 ++outputBits;
411 assert(outputBits <= 16);
412
413 parseVC5();
414 }
415
initVC5LogTable()416 void VC5Decompressor::initVC5LogTable() {
417 mVC5LogTable = decltype(mVC5LogTable)(
418 [outputBits = outputBits](unsigned i, unsigned tableSize) {
419 // The vanilla "inverse log" curve for decoding.
420 auto normalizedCurve = [](auto normalizedI) {
421 return (std::pow(113.0, normalizedI) - 1) / 112.0;
422 };
423
424 auto normalizeI = [tableSize](auto x) { return x / (tableSize - 1.0); };
425 auto denormalizeY = [maxVal = std::numeric_limits<uint16_t>::max()](
426 auto y) { return maxVal * y; };
427 // Adjust for output whitelevel bitdepth.
428 auto rescaleY = [outputBits](auto y) {
429 auto scale = 16 - outputBits;
430 return y >> scale;
431 };
432
433 const auto naiveY = denormalizeY(normalizedCurve(normalizeI(i)));
434 const auto intY = static_cast<unsigned int>(naiveY);
435 const auto rescaledY = rescaleY(intY);
436 return rescaledY;
437 });
438 }
439
parseVC5()440 void VC5Decompressor::parseVC5() {
441 mBs.setByteOrder(Endianness::big);
442
443 assert(mRaw->dim.x > 0);
444 assert(mRaw->dim.y > 0);
445
446 // All VC-5 data must start with "VC-%" (0x56432d35)
447 if (mBs.getU32() != 0x56432d35)
448 ThrowRDE("not a valid VC-5 datablock");
449
450 bool done = false;
451 while (!done) {
452 auto tag = static_cast<VC5Tag>(mBs.getU16());
453 uint16_t val = mBs.getU16();
454
455 bool optional = matches(tag, VC5Tag::Optional);
456 if (optional)
457 tag = -tag;
458
459 switch (tag) {
460 case VC5Tag::ChannelCount:
461 if (val != numChannels)
462 ThrowRDE("Bad channel count %u, expected %u", val, numChannels);
463 break;
464 case VC5Tag::ImageWidth:
465 if (val != mRaw->dim.x)
466 ThrowRDE("Image width mismatch: %u vs %u", val, mRaw->dim.x);
467 break;
468 case VC5Tag::ImageHeight:
469 if (val != mRaw->dim.y)
470 ThrowRDE("Image height mismatch: %u vs %u", val, mRaw->dim.y);
471 break;
472 case VC5Tag::LowpassPrecision:
473 if (val < PRECISION_MIN || val > PRECISION_MAX)
474 ThrowRDE("Invalid precision %i", val);
475 mVC5.lowpassPrecision = val;
476 break;
477 case VC5Tag::ChannelNumber:
478 if (val >= numChannels)
479 ThrowRDE("Bad channel number (%u)", val);
480 mVC5.iChannel = val;
481 break;
482 case VC5Tag::ImageFormat:
483 if (val != mVC5.imgFormat)
484 ThrowRDE("Image format %i is not 4(RAW)", val);
485 break;
486 case VC5Tag::SubbandCount:
487 if (val != numSubbands)
488 ThrowRDE("Unexpected subband count %u, expected %u", val, numSubbands);
489 break;
490 case VC5Tag::MaxBitsPerComponent:
491 if (val != VC5_LOG_TABLE_BITWIDTH) {
492 ThrowRDE("Bad bits per componend %u, not %u", val,
493 VC5_LOG_TABLE_BITWIDTH);
494 }
495 break;
496 case VC5Tag::PatternWidth:
497 if (val != mVC5.patternWidth)
498 ThrowRDE("Bad pattern width %u, not %u", val, mVC5.patternWidth);
499 break;
500 case VC5Tag::PatternHeight:
501 if (val != mVC5.patternHeight)
502 ThrowRDE("Bad pattern height %u, not %u", val, mVC5.patternHeight);
503 break;
504 case VC5Tag::SubbandNumber:
505 if (val >= numSubbands)
506 ThrowRDE("Bad subband number %u", val);
507 mVC5.iSubband = val;
508 break;
509 case VC5Tag::Quantization:
510 mVC5.quantization = static_cast<int16_t>(val);
511 break;
512 case VC5Tag::ComponentsPerSample:
513 if (val != mVC5.cps)
514 ThrowRDE("Bad component per sample count %u, not %u", val, mVC5.cps);
515 break;
516 case VC5Tag::PrescaleShift:
517 // FIXME: something is wrong. We get this before VC5Tag::ChannelNumber.
518 // Defaulting to 'mVC5.iChannel=0' seems to work *for existing samples*.
519 for (int iWavelet = 0; iWavelet < numWaveletLevels; ++iWavelet) {
520 auto& channel = channels[mVC5.iChannel];
521 auto& wavelet = channel.wavelets[1 + iWavelet];
522 wavelet.prescale =
523 extractHighBits(val, 2 * iWavelet, /*effectiveBitwidth=*/14) & 0x03;
524 }
525 break;
526 default: { // A chunk.
527 unsigned int chunkSize = 0;
528 if (matches(tag, VC5Tag::LARGE_CHUNK)) {
529 chunkSize = static_cast<unsigned int>(
530 ((static_cast<std::underlying_type<VC5Tag>::type>(tag) & 0xff)
531 << 16) |
532 (val & 0xffff));
533 } else if (matches(tag, VC5Tag::SMALL_CHUNK)) {
534 chunkSize = (val & 0xffff);
535 }
536
537 if (is(tag, VC5Tag::LargeCodeblock)) {
538 parseLargeCodeblock(mBs.getStream(chunkSize, 4));
539 break;
540 }
541
542 // And finally, we got here if we didn't handle this tag/maybe-chunk.
543
544 // Magic, all the other 'large' chunks are actually optional,
545 // and don't specify any chunk bytes-to-be-skipped.
546 if (matches(tag, VC5Tag::LARGE_CHUNK)) {
547 optional = true;
548 chunkSize = 0;
549 }
550
551 if (!optional) {
552 ThrowRDE("Unknown (unhandled) non-optional Tag 0x%04hx",
553 static_cast<std::underlying_type<VC5Tag>::type>(tag));
554 }
555
556 if (chunkSize)
557 mBs.skipBytes(chunkSize, 4);
558
559 break;
560 }
561 }
562
563 done = true;
564 for (int iChannel = 0; iChannel < numChannels && done; ++iChannel) {
565 Wavelet& wavelet = channels[iChannel].wavelets[0];
566 if (!wavelet.isBandValid(0))
567 done = false;
568 }
569 }
570 }
571
createDecodingTasks(ErrorLog & errLog,bool & exceptionThrown)572 void VC5Decompressor::Wavelet::AbstractDecodeableBand::createDecodingTasks(
573 ErrorLog& errLog, bool& exceptionThrown) noexcept {
574 auto& decodedData = data;
575
576 #ifdef HAVE_OPENMP
577 #pragma omp task default(none) shared(decodedData, errLog, exceptionThrown) \
578 depend(out \
579 : decodedData)
580 #endif
581 {
582 // Proceed only if decoding did not fail.
583 if (!readValue(exceptionThrown)) {
584 try {
585 assert(!decodedData.has_value() && "Decoded this band already?");
586 decodedData = decode();
587 } catch (RawspeedException& err) {
588 // Propagate the exception out of OpenMP magic.
589 errLog.setError(err.what());
590 #ifdef HAVE_OPENMP
591 #pragma omp atomic write
592 #endif
593 exceptionThrown = true;
594 }
595 }
596 }
597 }
598
LowPassBand(Wavelet & wavelet_,ByteStream bs_,uint16_t lowpassPrecision_)599 VC5Decompressor::Wavelet::LowPassBand::LowPassBand(Wavelet& wavelet_,
600 ByteStream bs_,
601 uint16_t lowpassPrecision_)
602 : AbstractDecodeableBand(wavelet_, std::move(bs_)),
603 lowpassPrecision(lowpassPrecision_) {
604 // Low-pass band is a uncompressed version of the image, hugely downscaled.
605 // It consists of width * height pixels, `lowpassPrecision` each.
606 // We can easily check that we have sufficient amount of bits to decode it.
607 const auto waveletArea = iPoint2D(wavelet.width, wavelet.height).area();
608 const auto bitsTotal = waveletArea * lowpassPrecision;
609 const auto bytesTotal = roundUpDivision(bitsTotal, 8);
610 bs = bs.getStream(bytesTotal); // And clamp the size while we are at it.
611 }
612
613 VC5Decompressor::BandData
decode() const614 VC5Decompressor::Wavelet::LowPassBand::decode() const noexcept {
615 BandData lowpass;
616 auto& band = lowpass.description;
617 band = Array2DRef<int16_t>::create(lowpass.storage, wavelet.width,
618 wavelet.height);
619
620 BitPumpMSB bits(bs);
621 for (auto row = 0; row < band.height; ++row) {
622 for (auto col = 0; col < band.width; ++col)
623 band(row, col) = static_cast<int16_t>(bits.getBits(lowpassPrecision));
624 }
625
626 return lowpass;
627 }
628
629 VC5Decompressor::BandData
decode() const630 VC5Decompressor::Wavelet::HighPassBand::decode() const {
631 class DeRLVer final {
632 BitPumpMSB bits;
633 const int16_t quant;
634
635 int16_t pixelValue = 0;
636 unsigned int numPixelsLeft = 0;
637
638 void decodeNextPixelGroup() {
639 assert(numPixelsLeft == 0);
640 std::tie(pixelValue, numPixelsLeft) = getRLV(bits);
641 }
642
643 public:
644 DeRLVer(const ByteStream& bs, int16_t quant_) : bits(bs), quant(quant_) {}
645
646 void verifyIsAtEnd() {
647 if (numPixelsLeft != 0)
648 ThrowRDE("Not all pixels consumed?");
649 decodeNextPixelGroup();
650 static_assert(decompand(MARKER_BAND_END) == MARKER_BAND_END,
651 "passthrough");
652 if (pixelValue != MARKER_BAND_END || numPixelsLeft != 0)
653 ThrowRDE("EndOfBand marker not found");
654 }
655
656 int16_t decode() {
657 auto dequantize = [quant = quant](int16_t val) -> int16_t {
658 return val * quant;
659 };
660
661 if (numPixelsLeft == 0) {
662 decodeNextPixelGroup();
663 pixelValue = dequantize(pixelValue);
664 }
665
666 if (numPixelsLeft == 0)
667 ThrowRDE("Got EndOfBand marker while looking for next pixel");
668
669 --numPixelsLeft;
670 return pixelValue;
671 }
672 };
673
674 // decode highpass band
675 DeRLVer d(bs, quant);
676 BandData highpass;
677 auto& band = highpass.description;
678 band = Array2DRef<int16_t>::create(highpass.storage, wavelet.width,
679 wavelet.height);
680 for (int row = 0; row != wavelet.height; ++row)
681 for (int col = 0; col != wavelet.width; ++col)
682 band(row, col) = d.decode();
683 d.verifyIsAtEnd();
684 return highpass;
685 }
686
parseLargeCodeblock(const ByteStream & bs)687 void VC5Decompressor::parseLargeCodeblock(const ByteStream& bs) {
688 static const auto subband_wavelet_index = []() {
689 std::array<int, numSubbands> wavelets;
690 int wavelet = 0;
691 for (auto i = wavelets.size() - 1; i > 0;) {
692 for (auto t = 0; t < numWaveletLevels; t++) {
693 wavelets[i] = wavelet;
694 i--;
695 }
696 if (i > 0)
697 wavelet++;
698 }
699 wavelets.front() = wavelet;
700 return wavelets;
701 }();
702 static const auto subband_band_index = []() {
703 std::array<int, numSubbands> bands;
704 bands.front() = 0;
705 for (auto i = 1U; i < bands.size();) {
706 for (int t = 1; t <= numWaveletLevels;) {
707 bands[i] = t;
708 t++;
709 i++;
710 }
711 }
712 return bands;
713 }();
714
715 if (!mVC5.iSubband.has_value())
716 ThrowRDE("Did not see VC5Tag::SubbandNumber yet");
717
718 const int idx = subband_wavelet_index[mVC5.iSubband.value()];
719 const int band = subband_band_index[mVC5.iSubband.value()];
720
721 auto& wavelets = channels[mVC5.iChannel].wavelets;
722
723 Wavelet& wavelet = wavelets[1 + idx];
724 if (wavelet.isBandValid(band)) {
725 ThrowRDE("Band %u for wavelet %u on channel %u was already seen", band, idx,
726 mVC5.iChannel);
727 }
728
729 std::unique_ptr<Wavelet::AbstractBand>& dstBand = wavelet.bands[band];
730 if (mVC5.iSubband.value() == 0) {
731 assert(band == 0);
732 // low-pass band, only one, for the smallest wavelet, per channel per image
733 if (!mVC5.lowpassPrecision.has_value())
734 ThrowRDE("Did not see VC5Tag::LowpassPrecision yet");
735 dstBand = std::make_unique<Wavelet::LowPassBand>(
736 wavelet, bs, mVC5.lowpassPrecision.value());
737 mVC5.lowpassPrecision.reset();
738 } else {
739 if (!mVC5.quantization.has_value())
740 ThrowRDE("Did not see VC5Tag::Quantization yet");
741 dstBand = std::make_unique<Wavelet::HighPassBand>(
742 wavelet, bs, mVC5.quantization.value());
743 mVC5.quantization.reset();
744 }
745 wavelet.setBandValid(band);
746
747 // If this wavelet is fully specified, mark the low-pass band of the
748 // next lower wavelet as specified.
749 if (wavelet.allBandsValid()) {
750 Wavelet& nextWavelet = wavelets[idx];
751 assert(!nextWavelet.isBandValid(0));
752 bool finalWavelet = idx == 0;
753 nextWavelet.bands[0] = std::make_unique<Wavelet::ReconstructableBand>(
754 wavelet, /*clampUint=*/finalWavelet, finalWavelet);
755 nextWavelet.setBandValid(0);
756 }
757
758 mVC5.iSubband.reset();
759 }
760
createWaveletBandDecodingTasks(bool & exceptionThrown) const761 void VC5Decompressor::createWaveletBandDecodingTasks(
762 bool& exceptionThrown) const noexcept {
763 for (int waveletLevel = numWaveletLevels; waveletLevel >= 0; waveletLevel--) {
764 const int numBandsInCurrentWavelet =
765 waveletLevel == 0 ? 1 : Wavelet::maxBands;
766 for (int bandId = 0; bandId != numBandsInCurrentWavelet; ++bandId) {
767 for (const auto& channel : channels) {
768 channel.wavelets[waveletLevel].bands[bandId]->createDecodingTasks(
769 static_cast<ErrorLog&>(*mRaw), exceptionThrown);
770 }
771 }
772 }
773 }
774
decodeThread(bool & exceptionThrown) const775 void VC5Decompressor::decodeThread(bool& exceptionThrown) const noexcept {
776 #ifdef HAVE_OPENMP
777 #pragma omp taskgroup
778 #pragma omp single
779 #endif
780 createWaveletBandDecodingTasks(exceptionThrown);
781
782 // Proceed only if decoding did not fail.
783 if (!readValue(exceptionThrown)) {
784 // And finally!
785 combineFinalLowpassBands();
786 }
787 }
788
decode(unsigned int offsetX,unsigned int offsetY,unsigned int width,unsigned int height)789 void VC5Decompressor::decode(unsigned int offsetX, unsigned int offsetY,
790 unsigned int width, unsigned int height) {
791 if (offsetX || offsetY || mRaw->dim != iPoint2D(width, height))
792 ThrowRDE("VC5Decompressor expects to fill the whole image, not some tile.");
793
794 initVC5LogTable();
795
796 alignas(RAWSPEED_CACHELINESIZE) bool exceptionThrown = false;
797
798 #ifdef HAVE_OPENMP
799 #pragma omp parallel default(none) shared(exceptionThrown) \
800 num_threads(rawspeed_get_number_of_processor_cores())
801 #endif
802 decodeThread(exceptionThrown);
803
804 std::string firstErr;
805 if (mRaw->isTooManyErrors(1, &firstErr)) {
806 assert(exceptionThrown);
807 ThrowRDE("Too many errors encountered. Giving up. First Error:\n%s",
808 firstErr.c_str());
809 } else {
810 assert(!exceptionThrown);
811 }
812 }
813
combineFinalLowpassBands() const814 void VC5Decompressor::combineFinalLowpassBands() const noexcept {
815 const Array2DRef<uint16_t> out(mRaw->getU16DataAsUncroppedArray2DRef());
816
817 const int width = out.width / 2;
818 const int height = out.height / 2;
819
820 assert(channels[0].wavelets[0].bands[0]->data.has_value() &&
821 channels[1].wavelets[0].bands[0]->data.has_value() &&
822 channels[2].wavelets[0].bands[0]->data.has_value() &&
823 channels[3].wavelets[0].bands[0]->data.has_value() &&
824 "Failed to reconstruct all final lowpass bands?");
825
826 const Array2DRef<const int16_t> lowbands0 =
827 channels[0].wavelets[0].bands[0]->data->description;
828 const Array2DRef<const int16_t> lowbands1 =
829 channels[1].wavelets[0].bands[0]->data->description;
830 const Array2DRef<const int16_t> lowbands2 =
831 channels[2].wavelets[0].bands[0]->data->description;
832 const Array2DRef<const int16_t> lowbands3 =
833 channels[3].wavelets[0].bands[0]->data->description;
834
835 // Convert to RGGB output
836 #ifdef HAVE_OPENMP
837 #pragma omp for schedule(static)
838 #endif
839 for (int row = 0; row < height; ++row) {
840 for (int col = 0; col < width; ++col) {
841 const int mid = 2048;
842
843 int gs = lowbands0(row, col);
844 int rg = lowbands1(row, col) - mid;
845 int bg = lowbands2(row, col) - mid;
846 int gd = lowbands3(row, col) - mid;
847
848 int r = gs + 2 * rg;
849 int b = gs + 2 * bg;
850 int g1 = gs + gd;
851 int g2 = gs - gd;
852
853 out(2 * row + 0, 2 * col + 0) = static_cast<uint16_t>(mVC5LogTable[r]);
854 out(2 * row + 0, 2 * col + 1) = static_cast<uint16_t>(mVC5LogTable[g1]);
855 out(2 * row + 1, 2 * col + 0) = static_cast<uint16_t>(mVC5LogTable[g2]);
856 out(2 * row + 1, 2 * col + 1) = static_cast<uint16_t>(mVC5LogTable[b]);
857 }
858 }
859 }
860
861 inline std::pair<int16_t /*value*/, unsigned int /*count*/>
getRLV(BitPumpMSB & bits)862 VC5Decompressor::getRLV(BitPumpMSB& bits) {
863 unsigned int iTab;
864
865 static constexpr auto maxBits = 1 + table17.entries[table17.length - 1].size;
866
867 // Ensure the maximum number of bits are cached to make peekBits() as fast as
868 // possible.
869 bits.fill(maxBits);
870 for (iTab = 0; iTab < table17.length; ++iTab) {
871 if (decompandedTable17[iTab].bits ==
872 bits.peekBitsNoFill(decompandedTable17[iTab].size))
873 break;
874 }
875 if (iTab >= table17.length)
876 ThrowRDE("Code not found in codebook");
877
878 bits.skipBitsNoFill(decompandedTable17[iTab].size);
879 int16_t value = decompandedTable17[iTab].value;
880 unsigned int count = decompandedTable17[iTab].count;
881 if (value != 0) {
882 if (bits.getBitsNoFill(1))
883 value = -(value);
884 }
885
886 return {value, count};
887 }
888
889 } // namespace rawspeed
890