View Javadoc
1   package nom.tam.fits.compression.algorithm.rice;
2   
3   import java.nio.Buffer;
4   import java.nio.ByteBuffer;
5   import java.nio.IntBuffer;
6   import java.nio.ShortBuffer;
7   import java.util.logging.Logger;
8   
9   import nom.tam.fits.compression.algorithm.api.ICompressor;
10  import nom.tam.fits.compression.algorithm.quant.QuantizeProcessor.DoubleQuantCompressor;
11  import nom.tam.fits.compression.algorithm.quant.QuantizeProcessor.FloatQuantCompressor;
12  import nom.tam.util.FitsIO;
13  import nom.tam.util.type.ElementType;
14  
15  /*
16   * #%L
17   * nom.tam FITS library
18   * %%
19   * Copyright (C) 1996 - 2024 nom-tam-fits
20   * %%
21   * This is free and unencumbered software released into the public domain.
22   *
23   * Anyone is free to copy, modify, publish, use, compile, sell, or
24   * distribute this software, either in source code form or as a compiled
25   * binary, for any purpose, commercial or non-commercial, and by any
26   * means.
27   *
28   * In jurisdictions that recognize copyright laws, the author or authors
29   * of this software dedicate any and all copyright interest in the
30   * software to the public domain. We make this dedication for the benefit
31   * of the public at large and to the detriment of our heirs and
32   * successors. We intend this dedication to be an overt act of
33   * relinquishment in perpetuity of all present and future rights to this
34   * software under copyright law.
35   *
36   * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
37   * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
38   * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
39   * IN NO EVENT SHALL THE AUTHORS BE LIABLE FOR ANY CLAIM, DAMAGES OR
40   * OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE,
41   * ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
42   * OTHER DEALINGS IN THE SOFTWARE.
43   * #L%
44   */
45  
46  /**
47   * (<i>for internal use</i>) The Rice compression algorithm. The original compression was designed by Rice, Yeh, and
48   * Miller the code was written by Richard White at STScI at the STScI and included (ported to c and adapted) in cfitsio
49   * by William Pence, NASA/GSFC. That code was then ported to java by R. van Nieuwenhoven. Later it was massively
50   * refactored to harmonize the different compression algorithms and reduce the duplicate code pieces without obscuring
51   * the algorithm itself as far as possible.
52   *
53   * @author     Richard White
54   * @author     William Pence
55   * @author     Richard van Nieuwenhoven
56   *
57   * @param  <T> the genetic type of NIO buffer on which this compressor operates.
58   */
59  @SuppressWarnings({"deprecation", "javadoc"})
60  public abstract class RiceCompressor<T extends Buffer> implements ICompressor<T> {
61  
62      public static class ByteRiceCompressor extends RiceCompressor<ByteBuffer> {
63  
64          private ByteBuffer pixelBuffer;
65  
66          /**
67           * Rice compression of byte streams with the default block size of 32.
68           * 
69           * @since  1.19.1
70           * 
71           * @author Attila Kovacs
72           */
73          public ByteRiceCompressor() {
74              this(new RiceCompressOption());
75          }
76  
77          public ByteRiceCompressor(RiceCompressOption option) {
78              super(option);
79          }
80  
81          @Override
82          public boolean compress(ByteBuffer buffer, ByteBuffer writeBuffer) {
83              pixelBuffer = buffer;
84              super.compress(buffer.limit(), pixelBuffer.get(pixelBuffer.position()), new BitBuffer(writeBuffer));
85              return true;
86          }
87  
88          @Override
89          public void decompress(ByteBuffer readBuffer, ByteBuffer buffer) {
90              pixelBuffer = buffer;
91              super.decompressBuffer(readBuffer, buffer.limit());
92          }
93  
94          @Override
95          protected int nextPixel() {
96              return pixelBuffer.get();
97          }
98  
99          @Override
100         protected void nextPixel(int pixel) {
101             pixelBuffer.put((byte) pixel);
102         }
103     }
104 
105     public static class DoubleRiceCompressor extends DoubleQuantCompressor {
106 
107         public DoubleRiceCompressor(RiceQuantizeCompressOption options) throws ClassCastException {
108             super(options, new IntRiceCompressor((RiceCompressOption) options.getCompressOption()));
109         }
110     }
111 
112     public static class FloatRiceCompressor extends FloatQuantCompressor {
113 
114         public FloatRiceCompressor(RiceQuantizeCompressOption options) throws ClassCastException {
115             super(options, new IntRiceCompressor((RiceCompressOption) options.getCompressOption()));
116         }
117     }
118 
119     public static class IntRiceCompressor extends RiceCompressor<IntBuffer> {
120 
121         private IntBuffer pixelBuffer;
122 
123         /**
124          * Rice compression of 32-bit integer streams with the default block size of 32.
125          * 
126          * @since  1.19.1
127          * 
128          * @author Attila Kovacs
129          */
130         public IntRiceCompressor() {
131             this(new RiceCompressOption());
132         }
133 
134         public IntRiceCompressor(RiceCompressOption option) {
135             super(option);
136         }
137 
138         @Override
139         public boolean compress(IntBuffer buffer, ByteBuffer writeBuffer) {
140             pixelBuffer = buffer;
141             super.compress(buffer.limit(), pixelBuffer.get(pixelBuffer.position()), new BitBuffer(writeBuffer));
142             return true;
143         }
144 
145         @Override
146         public void decompress(ByteBuffer readBuffer, IntBuffer buffer) {
147             pixelBuffer = buffer;
148             super.decompressBuffer(readBuffer, buffer.limit());
149         }
150 
151         @Override
152         protected int nextPixel() {
153             return pixelBuffer.get();
154         }
155 
156         @Override
157         protected void nextPixel(int pixel) {
158             pixelBuffer.put(pixel);
159         }
160     }
161 
162     public static class ShortRiceCompressor extends RiceCompressor<ShortBuffer> {
163 
164         private ShortBuffer pixelBuffer;
165 
166         /**
167          * Rice compression of 16-bit integer streams with the default block size of 32.
168          * 
169          * @since  1.19.1
170          * 
171          * @author Attila Kovacs
172          */
173         public ShortRiceCompressor() {
174             this(new RiceCompressOption());
175         }
176 
177         public ShortRiceCompressor(RiceCompressOption option) {
178             super(option);
179         }
180 
181         @Override
182         public boolean compress(ShortBuffer buffer, ByteBuffer writeBuffer) {
183             pixelBuffer = buffer;
184             super.compress(buffer.limit(), pixelBuffer.get(pixelBuffer.position()), new BitBuffer(writeBuffer));
185             return true;
186         }
187 
188         @Override
189         public void decompress(ByteBuffer readBuffer, ShortBuffer buffer) {
190             pixelBuffer = buffer;
191             super.decompressBuffer(readBuffer, buffer.limit());
192         }
193 
194         @Override
195         protected int nextPixel() {
196             return pixelBuffer.get();
197         }
198 
199         @Override
200         protected void nextPixel(int pixel) {
201             pixelBuffer.put((short) pixel);
202         }
203     }
204 
205     /**
206      * mask to convert a "unsigned" byte to a long.
207      */
208     private static final long UNSIGNED_BYTE_MASK = 0xFFL;
209 
210     /**
211      * mask to convert a "unsigned" short to a long.
212      */
213     private static final long UNSIGNED_SHORT_MASK = 0xFFFFL;
214 
215     /**
216      * mask to convert a "unsigned" int to a long.
217      */
218     private static final long UNSIGNED_INTEGER_MASK = 0xFFFFFFFFL;
219 
220     /**
221      * logger to log to.
222      */
223     private static final Logger LOG = Logger.getLogger(RiceCompressor.class.getName());
224 
225     private static final int BITS_OF_1_BYTE = 8;
226 
227     private static final int BITS_PER_BYTE = 8;
228 
229     private static final int BYTE_MASK = 0xff;
230 
231     private static final int FS_BITS_FOR_BYTE = 3;
232 
233     private static final int FS_BITS_FOR_INT = 5;
234 
235     private static final int FS_BITS_FOR_SHORT = 4;
236 
237     private static final int FS_MAX_FOR_BYTE = 6;
238 
239     private static final int FS_MAX_FOR_INT = 25;
240 
241     private static final int FS_MAX_FOR_SHORT = 14;
242 
243     /*
244      * nonzero_count is lookup table giving number of bits in 8-bit values not including leading zeros used in
245      * fits_rdecomp, fits_rdecomp_short and fits_rdecomp_byte.
246      *
247      * @formatter:off
248      */
249     private static final int[] NONZERO_COUNT = {0, 1, 2, 2, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 5, 5, 5, 5,
250             5, 5, 5, 5, 5, 5, 5, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6,
251             6, 6, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7,
252             7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 8, 8, 8, 8, 8, 8, 8, 8,
253             8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8,
254             8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8,
255             8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8,
256             8, 8, 8, 8, 8, 8, 8, 8, 8};
257     // @formatter:on
258 
259     private final int bBits;
260 
261     private final int bitsPerPixel;
262 
263     private final int blockSize;
264 
265     private final int fsBits;
266 
267     private final int fsMax;
268 
269     private RiceCompressor(RiceCompressOption option) throws UnsupportedOperationException {
270         blockSize = option.getBlockSize();
271         if (option.getBytePix() == ElementType.BYTE.size()) {
272             fsBits = FS_BITS_FOR_BYTE;
273             fsMax = FS_MAX_FOR_BYTE;
274             bitsPerPixel = FitsIO.BITS_OF_1_BYTE;
275         } else if (option.getBytePix() == ElementType.SHORT.size()) {
276             fsBits = FS_BITS_FOR_SHORT;
277             fsMax = FS_MAX_FOR_SHORT;
278             bitsPerPixel = FitsIO.BITS_OF_2_BYTES;
279         } else if (option.getBytePix() == ElementType.INT.size()) {
280             fsBits = FS_BITS_FOR_INT;
281             fsMax = FS_MAX_FOR_INT;
282             bitsPerPixel = FitsIO.BITS_OF_4_BYTES;
283         } else {
284             throw new UnsupportedOperationException("Implemented for 1/2/4 bytes only");
285         }
286         /*
287          * From bsize derive: FSBITS = # bits required to store FS FSMAX = maximum value for FS BBITS = bits/pixel for
288          * direct coding
289          */
290         bBits = 1 << fsBits;
291     }
292 
293     /**
294      * <p>
295      * undo mapping and differencing Note that some of these operations will overflow the unsigned int arithmetic --
296      * that's OK, it all works out to give the right answers in the output file.
297      * </p>
298      * <p>
299      * In java this is more complicated because of the missing unsigned integers. trying to simulate the behavior
300      * </p>
301      *
302      * @param  lastpix the current last pix value
303      * @param  diff    the difference to "add"
304      *
305      * @return         return the new lastpiy value
306      */
307     private long undoMappingAndDifferencing(long lastpix, long diff) {
308         diff &= UNSIGNED_INTEGER_MASK;
309         if ((diff & 1) == 0) {
310             diff = diff >>> 1;
311         } else {
312             diff = diff >>> 1 ^ UNSIGNED_INTEGER_MASK;
313         }
314         lastpix = diff + lastpix & UNSIGNED_INTEGER_MASK;
315         nextPixel((int) lastpix);
316         return lastpix;
317     }
318 
319     /**
320      * compress the integer tiledImageOperation on a rise compressed byte buffer.
321      *
322      * @param dataLength length of the data to compress
323      * @param firstPixel the value of the first pixel
324      * @param buffer     the buffer to write to
325      */
326     protected void compress(final int dataLength, int firstPixel, BitBuffer buffer) {
327         /* the first difference will always be zero */
328         int lastpix = firstPixel;
329         /* write out first int value to the first 4 bytes of the buffer */
330         buffer.putInt(firstPixel, bitsPerPixel);
331         int thisblock = blockSize;
332         for (int i = 0; i < dataLength; i += blockSize) {
333             /* last block may be shorter */
334             if (dataLength - i < blockSize) {
335                 thisblock = dataLength - i;
336             }
337             /*
338              * Compute differences of adjacent pixels and map them to unsigned values. Note that this may overflow the
339              * integer variables -- that's OK, because we can recover when decompressing. If we were compressing shorts
340              * or bytes, would want to do this arithmetic with short/byte working variables (though diff will still be
341              * passed as an int.) compute sum of mapped pixel values at same time use double precision for sum to allow
342              * 32-bit integer inputs
343              */
344             long[] diff = new long[blockSize];
345             double pixelsum = 0.0;
346             int nextpix;
347             /*
348              * tiledImageOperation for differences mapped to non-negative values
349              */
350             for (int j = 0; j < thisblock; j++) {
351                 nextpix = nextPixel();
352                 long pdiff = (nextpix - lastpix);
353                 diff[j] = (pdiff < 0 ? (pdiff << 1) ^ UNSIGNED_INTEGER_MASK : pdiff << 1) & UNSIGNED_INTEGER_MASK;
354                 pixelsum += diff[j];
355                 lastpix = nextpix;
356             }
357 
358             /*
359              * compute number of bits to split from sum
360              */
361             double dpsum = (pixelsum - thisblock / 2d - 1d) / thisblock;
362             if (dpsum < 0) {
363                 dpsum = 0.0;
364             }
365             long psum = (long) dpsum >> 1;
366             int fs;
367             for (fs = 0; psum > 0; fs++) { // NOSONAR
368                 psum >>= 1;
369             }
370 
371             /*
372              * write the codes fsbits ID bits used to indicate split level
373              */
374             if (fs >= fsMax) {
375                 /*
376                  * Special high entropy case when FS >= fsmax Just write pixel difference values directly, no Rice
377                  * coding at all.
378                  */
379                 buffer.putInt(fsMax + 1, fsBits);
380                 for (int j = 0; j < thisblock; j++) {
381                     buffer.putLong(diff[j], bBits);
382                 }
383             } else if (fs == 0 && pixelsum == 0) { // NOSONAR
384                 /*
385                  * special low entropy case when FS = 0 and pixelsum=0 (all pixels in block are zero.) Output a 0 and
386                  * return
387                  */
388                 buffer.putInt(0, fsBits);
389             } else {
390                 /* normal case: not either very high or very low entropy */
391                 buffer.putInt(fs + 1, fsBits);
392                 int fsmask = (1 << fs) - 1;
393                 /*
394                  * local copies of bit buffer to improve optimization
395                  */
396                 int bitsToGo = buffer.missingBitsInCurrentByte();
397                 int bitBuffer = buffer.bitbuffer() >> bitsToGo;
398                 buffer.movePosition(bitsToGo - BITS_OF_1_BYTE);
399                 for (int j = 0; j < thisblock; j++) {
400                     int v = (int) diff[j];
401                     int top = v >> fs;
402                     /*
403                      * top is coded by top zeros + 1
404                      */
405                     if (bitsToGo >= top + 1) {
406                         bitBuffer <<= top + 1;
407                         bitBuffer |= 1;
408                         bitsToGo -= top + 1;
409                     } else {
410                         bitBuffer <<= bitsToGo;
411                         buffer.putByte((byte) (bitBuffer & BYTE_MASK));
412                         for (top -= bitsToGo; top >= BITS_OF_1_BYTE; top -= BITS_OF_1_BYTE) {
413                             buffer.putByte((byte) 0);
414                         }
415                         bitBuffer = 1;
416                         bitsToGo = BITS_OF_1_BYTE - 1 - top;
417                     }
418                     /*
419                      * bottom FS bits are written without coding code is output_nbits, moved into this routine to reduce
420                      * overheads This code potentially breaks if FS>24, so I am limiting FS to 24 by choice of FSMAX
421                      * above.
422                      */
423                     if (fs > 0) {
424                         bitBuffer <<= fs;
425                         bitBuffer |= v & fsmask;
426                         bitsToGo -= fs;
427                         while (bitsToGo <= 0) {
428                             buffer.putByte((byte) (bitBuffer >> -bitsToGo & BYTE_MASK));
429                             bitsToGo += BITS_OF_1_BYTE;
430                         }
431                     }
432                 }
433                 buffer.putByte((byte) (bitBuffer & BYTE_MASK), BITS_OF_1_BYTE - bitsToGo);
434             }
435         }
436         buffer.close();
437     }
438 
439     /**
440      * decompress the readbuffer and fill the pixelarray.
441      *
442      * @param readBuffer input buffer
443      * @param nx         the number of pixel to uncompress
444      */
445     protected void decompressBuffer(final ByteBuffer readBuffer, final int nx) {
446         /* first x bytes of input buffer contain the value of the first */
447         /* x byte integer value, without any encoding */
448         long lastpix = 0L;
449         if (bitsPerPixel == ElementType.BYTE.bitPix()) {
450             lastpix = readBuffer.get() & UNSIGNED_BYTE_MASK;
451         } else if (bitsPerPixel == ElementType.SHORT.bitPix()) {
452             lastpix = readBuffer.getShort() & UNSIGNED_SHORT_MASK;
453         } else {
454             // Must be (this.bitsPerPixel == ElementType.INT.bitPix())
455             lastpix = readBuffer.getInt() & UNSIGNED_INTEGER_MASK;
456         }
457         long b = readBuffer.get() & BYTE_MASK; /* bit buffer */
458         int nbits = BITS_PER_BYTE; /* number of bits remaining in b */
459         for (int i = 0; i < nx;) {
460             /* get the FS value from first fsbits */
461             nbits -= fsBits;
462             while (nbits < 0) {
463                 b = b << BITS_PER_BYTE | readBuffer.get() & BYTE_MASK;
464                 nbits += BITS_PER_BYTE;
465             }
466             long fs = (b >>> nbits) - 1L;
467 
468             b &= (1 << nbits) - 1;
469             /* loop over the next block */
470             int imax = i + blockSize;
471             if (imax > nx) {
472                 imax = nx;
473             }
474             if (fs < 0) {
475                 /* low-entropy case, all zero differences */
476                 for (; i < imax; i++) {
477                     nextPixel((int) lastpix);
478                 }
479             } else if (fs == fsMax) {
480                 /* high-entropy case, directly coded pixel values */
481                 for (; i < imax; i++) {
482                     int k = bBits - nbits;
483                     long diff = b << k;
484                     for (k -= BITS_PER_BYTE; k >= 0; k -= BITS_PER_BYTE) {
485                         b = readBuffer.get() & BYTE_MASK;
486                         diff |= b << k;
487                     }
488                     if (nbits > 0) {
489                         b = readBuffer.get() & BYTE_MASK;
490                         diff |= b >>> -k;
491                         b &= (1 << nbits) - 1L;
492                     } else {
493                         b = 0;
494                     }
495                     lastpix = undoMappingAndDifferencing(lastpix, diff);
496                 }
497             } else {
498                 /* normal case, Rice coding */
499                 for (; i < imax; i++) {
500                     /* count number of leading zeros */
501                     while (b == 0) {
502                         nbits += BITS_PER_BYTE;
503                         b = readBuffer.get() & BYTE_MASK;
504                     }
505                     long nzero = nbits - NONZERO_COUNT[(int) (b & BYTE_MASK)];
506                     nbits -= nzero + 1;
507                     /* flip the leading one-bit */
508                     b ^= 1 << nbits;
509                     /* get the FS trailing bits */
510                     nbits -= fs;
511                     while (nbits < 0) {
512                         b = b << BITS_PER_BYTE | readBuffer.get() & BYTE_MASK;
513                         nbits += BITS_PER_BYTE;
514                     }
515                     long diff = nzero << fs | b >> nbits;
516                     b &= (1 << nbits) - 1L;
517 
518                     lastpix = undoMappingAndDifferencing(lastpix, diff);
519                 }
520             }
521         }
522         if (readBuffer.limit() > readBuffer.position()) {
523             LOG.warning("decompressing left over some extra bytes got: " + readBuffer.limit() + " but needed only "
524                     + readBuffer.position());
525         }
526 
527     }
528 
529     protected abstract int nextPixel();
530 
531     protected abstract void nextPixel(int pixel);
532 
533 }