1 /*
2 RawSpeed - RAW file decoder.
3
4 Copyright (C) 2017 Vasily Khoruzhick
5
6 This library is free software; you can redistribute it and/or
7 modify it under the terms of the GNU Lesser General Public
8 License as published by the Free Software Foundation; either
9 version 2 of the License, or (at your option) any later version.
10
11 This library is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 Lesser General Public License for more details.
15
16 You should have received a copy of the GNU Lesser General Public
17 License along with this library; if not, write to the Free Software
18 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
19 */
20
21 #include "rawspeedconfig.h"
22
23 #ifdef HAVE_ZLIB
24
25 #include "common/Point.h" // for iPoint2D
26 #include "decoders/RawDecoderException.h" // for ThrowRDE
27 #include "decompressors/DeflateDecompressor.h"
28 #include "io/Endianness.h" // for getHostEndianness, Endianness
29 #include <cassert> // for assert
30 #include <cstdint> // for uint32_t, uint16_t
31 #include <cstdio> // for size_t
32 #include <zlib.h> // for uncompress, zError, Z_OK
33
34 namespace rawspeed {
35
36 // decodeFPDeltaRow(): MIT License, copyright 2014 Javier Celaya
37 // <jcelaya@gmail.com>
decodeFPDeltaRow(unsigned char * src,unsigned char * dst,size_t tileWidth,size_t realTileWidth,unsigned int bytesps,int factor)38 static inline void decodeFPDeltaRow(unsigned char* src, unsigned char* dst,
39 size_t tileWidth, size_t realTileWidth,
40 unsigned int bytesps, int factor) {
41 // DecodeDeltaBytes
42 for (size_t col = factor; col < realTileWidth * bytesps; ++col) {
43 // Yes, this is correct, and is symmetrical with EncodeDeltaBytes in
44 // hdrmerge, and they both combined are lossless.
45 // This is indeed working in modulo-2^n arighmetics.
46 src[col] = static_cast<unsigned char>(src[col] + src[col - factor]);
47 }
48 // Reorder bytes into the image
49 // 16 and 32-bit versions depend on local architecture, 24-bit does not
50 if (bytesps == 3) {
51 for (size_t col = 0; col < tileWidth; ++col) {
52 dst[col * 3] = src[col];
53 dst[col * 3 + 1] = src[col + realTileWidth];
54 dst[col * 3 + 2] = src[col + realTileWidth * 2];
55 }
56 } else {
57 if (getHostEndianness() == Endianness::little) {
58 for (size_t col = 0; col < tileWidth; ++col) {
59 for (size_t byte = 0; byte < bytesps; ++byte)
60 dst[col * bytesps + byte] =
61 src[col + realTileWidth * (bytesps - byte - 1)];
62 }
63 } else {
64 for (size_t col = 0; col < tileWidth; ++col) {
65 for (size_t byte = 0; byte < bytesps; ++byte)
66 dst[col * bytesps + byte] = src[col + realTileWidth * byte];
67 }
68 }
69 }
70 }
71
fp16ToFloat(uint16_t fp16)72 static inline uint32_t __attribute__((const)) fp16ToFloat(uint16_t fp16) {
73 // IEEE-754-2008: binary16:
74 // bit 15 - sign
75 // bits 14-10 - exponent (5 bit)
76 // bits 9-0 - fraction (10 bit)
77 //
78 // exp = 0, fract = +-0: zero
79 // exp = 0; fract != 0: subnormal numbers
80 // equation: -1 ^ sign * 2 ^ -14 * 0.fraction
81 // exp = 1..30: normalized value
82 // equation: -1 ^ sign * 2 ^ (exponent - 15) * 1.fraction
83 // exp = 31, fract = +-0: +-infinity
84 // exp = 31, fract != 0: NaN
85
86 uint32_t sign = (fp16 >> 15) & 1;
87 uint32_t fp16_exponent = (fp16 >> 10) & ((1 << 5) - 1);
88 uint32_t fp16_fraction = fp16 & ((1 << 10) - 1);
89
90 // Normalized or zero
91 // binary32 equation: -1 ^ sign * 2 ^ (exponent - 127) * 1.fraction
92 // => exponent32 - 127 = exponent16 - 15, exponent32 = exponent16 + 127 - 15
93 uint32_t fp32_exponent = fp16_exponent + 127 - 15;
94 uint32_t fp32_fraction = fp16_fraction
95 << (23 - 10); // 23 is binary32 fraction size
96
97 if (fp16_exponent == 31) {
98 // Infinity or NaN
99 fp32_exponent = 255;
100 } else if (fp16_exponent == 0) {
101 if (fp16_fraction == 0) {
102 // +-Zero
103 fp32_exponent = 0;
104 fp32_fraction = 0;
105 } else {
106 // Subnormal numbers
107 // binary32 equation: -1 ^ sign * 2 ^ (exponent - 127) * 1.fraction
108 // binary16 equation: -1 ^ sign * 2 ^ -14 * 0.fraction, we can represent
109 // it as a normalized value in binary32, we have to shift fraction until
110 // we get 1.new_fraction and decrement exponent for each shift
111 fp32_exponent = -14 + 127;
112 while (!(fp32_fraction & (1 << 23))) {
113 fp32_exponent -= 1;
114 fp32_fraction <<= 1;
115 }
116 fp32_fraction &= ((1 << 23) - 1);
117 }
118 }
119 return (sign << 31) | (fp32_exponent << 23) | fp32_fraction;
120 }
121
fp24ToFloat(uint32_t fp24)122 static inline uint32_t __attribute__((const)) fp24ToFloat(uint32_t fp24) {
123 // binary24: Not a part of IEEE754-2008, but format is obvious,
124 // see https://en.wikipedia.org/wiki/Minifloat
125 // bit 23 - sign
126 // bits 22-16 - exponent (7 bit)
127 // bits 15-0 - fraction (16 bit)
128 //
129 // exp = 0, fract = +-0: zero
130 // exp = 0; fract != 0: subnormal numbers
131 // equation: -1 ^ sign * 2 ^ -62 * 0.fraction
132 // exp = 1..126: normalized value
133 // equation: -1 ^ sign * 2 ^ (exponent - 63) * 1.fraction
134 // exp = 127, fract = +-0: +-infinity
135 // exp = 127, fract != 0: NaN
136
137 uint32_t sign = (fp24 >> 23) & 1;
138 uint32_t fp24_exponent = (fp24 >> 16) & ((1 << 7) - 1);
139 uint32_t fp24_fraction = fp24 & ((1 << 16) - 1);
140
141 // Normalized or zero
142 // binary32 equation: -1 ^ sign * 2 ^ (exponent - 127) * 1.fraction
143 // => exponent32 - 127 = exponent24 - 64, exponent32 = exponent16 + 127 - 63
144 uint32_t fp32_exponent = fp24_exponent + 127 - 63;
145 uint32_t fp32_fraction = fp24_fraction
146 << (23 - 16); // 23 is binary 32 fraction size
147
148 if (fp24_exponent == 127) {
149 // Infinity or NaN
150 fp32_exponent = 255;
151 } else if (fp24_exponent == 0) {
152 if (fp24_fraction == 0) {
153 // +-Zero
154 fp32_exponent = 0;
155 fp32_fraction = 0;
156 } else {
157 // Subnormal numbers
158 // binary32 equation: -1 ^ sign * 2 ^ (exponent - 127) * 1.fraction
159 // binary24 equation: -1 ^ sign * 2 ^ -62 * 0.fraction, we can represent
160 // it as a normalized value in binary32, we have to shift fraction until
161 // we get 1.new_fraction and decrement exponent for each shift
162 fp32_exponent = -62 + 127;
163 while (!(fp32_fraction & (1 << 23))) {
164 fp32_exponent -= 1;
165 fp32_fraction <<= 1;
166 }
167 fp32_fraction &= ((1 << 23) - 1);
168 }
169 }
170 return (sign << 31) | (fp32_exponent << 23) | fp32_fraction;
171 }
172
expandFP16(unsigned char * dst,int width)173 static inline void expandFP16(unsigned char* dst, int width) {
174 auto* dst16 = reinterpret_cast<uint16_t*>(dst);
175 auto* dst32 = reinterpret_cast<uint32_t*>(dst);
176
177 for (int x = width - 1; x >= 0; x--)
178 dst32[x] = fp16ToFloat(dst16[x]);
179 }
180
expandFP24(unsigned char * dst,int width)181 static inline void expandFP24(unsigned char* dst, int width) {
182 auto* dst32 = reinterpret_cast<uint32_t*>(dst);
183 dst += (width - 1) * 3;
184 for (int x = width - 1; x >= 0; x--) {
185 dst32[x] = fp24ToFloat((dst[0] << 16) | (dst[1] << 8) | dst[2]);
186 dst -= 3;
187 }
188 }
189
decode(std::unique_ptr<unsigned char[]> * uBuffer,iPoint2D maxDim,iPoint2D dim,iPoint2D off)190 void DeflateDecompressor::decode(
191 std::unique_ptr<unsigned char[]>* uBuffer, // NOLINT
192 iPoint2D maxDim, iPoint2D dim, iPoint2D off) {
193 uLongf dstLen = sizeof(float) * maxDim.area();
194
195 if (!*uBuffer)
196 *uBuffer =
197 std::unique_ptr<unsigned char[]>(new unsigned char[dstLen]); // NOLINT
198
199 const auto cSize = input.getRemainSize();
200 const unsigned char* cBuffer = input.getData(cSize);
201
202 int err = uncompress(uBuffer->get(), &dstLen, cBuffer, cSize);
203 if (err != Z_OK) {
204 ThrowRDE("failed to uncompress tile: %d (%s)", err, zError(err));
205 }
206
207 int predFactor = 0;
208 switch (predictor) {
209 case 3:
210 predFactor = 1;
211 break;
212 case 34894:
213 predFactor = 2;
214 break;
215 case 34895:
216 predFactor = 4;
217 break;
218 default:
219 predFactor = 0;
220 break;
221 }
222 predFactor *= mRaw->getCpp();
223
224 int bytesps = bps / 8;
225
226 for (auto row = 0; row < dim.y; ++row) {
227 unsigned char* src = uBuffer->get() + row * maxDim.x * bytesps;
228 unsigned char* dst = static_cast<unsigned char*>(mRaw->getData()) +
229 ((off.y + row) * mRaw->pitch + off.x * sizeof(float));
230
231 if (predFactor)
232 decodeFPDeltaRow(src, dst, dim.x, maxDim.x, bytesps, predFactor);
233
234 assert(bytesps >= 2 && bytesps <= 4);
235 switch (bytesps) {
236 case 2:
237 expandFP16(dst, dim.x);
238 break;
239 case 3:
240 expandFP24(dst, dim.x);
241 break;
242 case 4:
243 // No need to expand FP32
244 break;
245 default:
246 __builtin_unreachable();
247 }
248 }
249 }
250
251 } // namespace rawspeed
252
253 #else
254
255 #pragma message \
256 "ZLIB is not present! Deflate compression will not be supported!"
257
258 #endif
259