001/*
002 * $RCSfile: StdQuantizer.java,v $
003 * $Revision: 1.1 $
004 * $Date: 2005/02/11 05:02:20 $
005 * $State: Exp $
006 *
007 * Class:                   StdQuantizer
008 *
009 * Description:             Scalar deadzone quantizer of integer or float
010 *                          data.
011 *
012 *                          Mergerd from StdQuantizerInt and
013 *                          StdQuantizerFloat from Joel Askelof.
014 *
015 *
016 * COPYRIGHT:
017 *
018 * This software module was originally developed by Raphaël Grosbois and
019 * Diego Santa Cruz (Swiss Federal Institute of Technology-EPFL); Joel
020 * Askelöf (Ericsson Radio Systems AB); and Bertrand Berthelot, David
021 * Bouchard, Félix Henry, Gerard Mozelle and Patrice Onno (Canon Research
022 * Centre France S.A) in the course of development of the JPEG2000
023 * standard as specified by ISO/IEC 15444 (JPEG 2000 Standard). This
024 * software module is an implementation of a part of the JPEG 2000
025 * Standard. Swiss Federal Institute of Technology-EPFL, Ericsson Radio
026 * Systems AB and Canon Research Centre France S.A (collectively JJ2000
027 * Partners) agree not to assert against ISO/IEC and users of the JPEG
028 * 2000 Standard (Users) any of their rights under the copyright, not
029 * including other intellectual property rights, for this software module
030 * with respect to the usage by ISO/IEC and Users of this software module
031 * or modifications thereof for use in hardware or software products
032 * claiming conformance to the JPEG 2000 Standard. Those intending to use
033 * this software module in hardware or software products are advised that
034 * their use may infringe existing patents. The original developers of
035 * this software module, JJ2000 Partners and ISO/IEC assume no liability
036 * for use of this software module or modifications thereof. No license
037 * or right to this software module is granted for non JPEG 2000 Standard
038 * conforming products. JJ2000 Partners have full right to use this
039 * software module for his/her own purpose, assign or donate this
040 * software module to any third party and to inhibit third parties from
041 * using this software module for non JPEG 2000 Standard conforming
042 * products. This copyright notice must be included in all copies or
043 * derivative works of this software module.
044 *
045 * Copyright (c) 1999/2000 JJ2000 Partners.
046 * */
047package jj2000.j2k.quantization.quantizer;
048import jj2000.j2k.codestream.writer.*;
049import jj2000.j2k.wavelet.analysis.*;
050import jj2000.j2k.quantization.*;
051import jj2000.j2k.wavelet.*;
052import jj2000.j2k.image.*;
053import jj2000.j2k.*;
054
055import com.sun.media.imageioimpl.plugins.jpeg2000.J2KImageWriteParamJava;
056
057/**
058 * This class implements scalar quantization of integer or floating-point
059 * valued source data. The source data is the wavelet transformed image data
060 * and the output is the quantized wavelet coefficients represented in
061 * sign-magnitude (see below).
062 *
063 * <P>Sign magnitude representation is used (instead of two's complement) for
064 * the output data. The most significant bit is used for the sign (0 if
065 * positive, 1 if negative). Then the magnitude of the quantized coefficient
066 * is stored in the next M most significat bits. The rest of the bits (least
067 * significant bits) can contain a fractional value of the quantized
068 * coefficient. This fractional value is not to be coded by the entropy
069 * coder. However, it can be used to compute rate-distortion measures with
070 * greater precision.
071 *
072 * <P>The value of M is determined for each subband as the sum of the number
073 * of guard bits G and the nominal range of quantized wavelet coefficients in
074 * the corresponding subband (Rq), minus 1:
075 *
076 * <P>M = G + Rq -1
077 *
078 * <P>The value of G should be the same for all subbands. The value of Rq
079 * depends on the quantization step size, the nominal range of the component
080 * before the wavelet transform and the analysis gain of the subband (see
081 * Subband).
082 *
083 * <P>The blocks of data that are requested should not cross subband
084 * boundaries.
085 *
086 * @see Subband
087 *
088 * @see Quantizer
089 * */
090public class StdQuantizer extends Quantizer {
091
092    /** The number of mantissa bits for the quantization steps */
093    public final static int QSTEP_MANTISSA_BITS = 11;
094
095    /** The number of exponent bits for the quantization steps */
096    // NOTE: formulas in 'convertFromExpMantissa()' and
097    // 'convertToExpMantissa()' methods do not support more than 5 bits.
098    public final static int QSTEP_EXPONENT_BITS = 5;
099
100    /** The maximum value of the mantissa for the quantization steps */
101    public final static int QSTEP_MAX_MANTISSA = (1<<QSTEP_MANTISSA_BITS)-1;
102
103    /** The maximum value of the exponent for the quantization steps */
104    public final static int QSTEP_MAX_EXPONENT = (1<<QSTEP_EXPONENT_BITS)-1;
105
106    /** Natural log of 2, used as a convenience variable */
107    private static double log2 = Math.log(2);
108
109    /** The quantization type specifications */
110    private QuantTypeSpec qts;
111
112    /** The quantization step size specifications */
113    private QuantStepSizeSpec qsss;
114
115    /** The guard bits specifications */
116    private GuardBitsSpec gbs;
117
118    /** The 'CBlkWTDataFloat' object used to request data, used when
119     * quantizing floating-point data. */
120    // This variable makes the class thread unsafe, but it avoids allocating
121    // new objects for code-block that is quantized.
122    private CBlkWTDataFloat infblk;
123
124    /**
125     * Initializes the source of wavelet transform coefficients. The
126     * constructor takes information on whether the quantizer is in
127     * reversible, derived or expounded mode. If the quantizer is reversible
128     * the value of 'derived' is ignored. If the source data is not integer
129     * (int) then the quantizer can not be reversible.
130     *
131     * <P> After initializing member attributes, getAnSubbandTree is called for
132     * all components setting the 'stepWMSE' for all subbands in the current
133     * tile.
134     *
135     * @param src The source of wavelet transform coefficients.
136     *
137     * @param encSpec The encoder specifications
138     * */
139    public StdQuantizer(CBlkWTDataSrc src,J2KImageWriteParamJava wp){
140        super(src);
141        qts  = wp.getQuantizationType();
142        qsss = wp.getQuantizationStep();
143        gbs  = wp.getGuardBits();
144    }
145
146    /**
147     * Returns the quantization type spec object associated to the quantizer.
148     *
149     * @return The quantization type spec
150     * */
151    public QuantTypeSpec getQuantTypeSpec(){
152        return qts;
153    }
154
155    /**
156     * Returns the number of guard bits used by this quantizer in the given
157     * tile-component.
158     *
159     * @param t Tile index
160     *
161     * @param c Component index
162     *
163     * @return The number of guard bits
164     * */
165    public int getNumGuardBits(int t,int c){
166        return ((Integer)gbs.getTileCompVal(t,c)).intValue();
167    }
168
169    /**
170     * Returns true if the quantized data is reversible, for the specified
171     * tile-component. For the quantized data to be reversible it is necessary
172     * and sufficient that the quantization is reversible.
173     *
174     * @param t The tile to test for reversibility
175     *
176     * @param c The component to test for reversibility
177     *
178     * @return True if the quantized data is reversible, false if not.
179     * */
180    public boolean isReversible(int t,int c){
181        return qts.isReversible(t,c);
182    }
183
184    /**
185     * Returns true if given tile-component uses derived quantization step
186     * sizes.
187     *
188     * @param t Tile index
189     *
190     * @param c Component index
191     *
192     * @return True if derived
193     *
194     */
195    public boolean isDerived(int t,int c){
196        return qts.isDerived(t,c);
197    }
198
199    /**
200     * Returns the next code-block in the current tile for the specified
201     * component, as a copy (see below). The order in which code-blocks are
202     * returned is not specified. However each code-block is returned only
203     * once and all code-blocks will be returned if the method is called 'N'
204     * times, where 'N' is the number of code-blocks in the tile. After all
205     * the code-blocks have been returned for the current tile calls to this
206     * method will return 'null'.
207     *
208     * <P>When changing the current tile (through 'setTile()' or 'nextTile()')
209     * this method will always return the first code-block, as if this method
210     * was never called before for the new current tile.
211     *
212     * <P>The data returned by this method is always a copy of the
213     * data. Therfore it can be modified "in place" without any problems after
214     * being returned. The 'offset' of the returned data is 0, and the 'scanw'
215     * is the same as the code-block width. See the 'CBlkWTData' class.
216     *
217     * <P>The 'ulx' and 'uly' members of the returned 'CBlkWTData' object
218     * contain the coordinates of the top-left corner of the block, with
219     * respect to the tile, not the subband.
220     *
221     * @param c The component for which to return the next code-block.
222     *
223     * @param cblk If non-null this object will be used to return the new
224     * code-block. If null a new one will be allocated and returned. If the
225     * "data" array of the object is non-null it will be reused, if possible,
226     * to return the data.
227     *
228     * @return The next code-block in the current tile for component 'n', or
229     * null if all code-blocks for the current tile have been returned.
230     *
231     * @see CBlkWTData
232     * */
233    public CBlkWTData getNextCodeBlock(int c,CBlkWTData cblk) {
234        return getNextInternCodeBlock(c,cblk);
235    }
236
237    /**
238     * Returns the next code-block in the current tile for the specified
239     * component. The order in which code-blocks are returned is not
240     * specified. However each code-block is returned only once and all
241     * code-blocks will be returned if the method is called 'N' times, where
242     * 'N' is the number of code-blocks in the tile. After all the code-blocks
243     * have been returned for the current tile calls to this method will
244     * return 'null'.
245     *
246     * <P>When changing the current tile (through 'setTile()' or 'nextTile()')
247     * this method will always return the first code-block, as if this method
248     * was never called before for the new current tile.
249     *
250     * <P>The data returned by this method can be the data in the internal
251     * buffer of this object, if any, and thus can not be modified by the
252     * caller. The 'offset' and 'scanw' of the returned data can be
253     * arbitrary. See the 'CBlkWTData' class.
254     *
255     * <P>The 'ulx' and 'uly' members of the returned 'CBlkWTData' object
256     * contain the coordinates of the top-left corner of the block, with
257     * respect to the tile, not the subband.
258     *
259     * @param c The component for which to return the next code-block.
260     *
261     * @param cblk If non-null this object will be used to return the new
262     * code-block. If null a new one will be allocated and returned. If the
263     * "data" array of the object is non-null it will be reused, if possible,
264     * to return the data.
265     *
266     * @return The next code-block in the current tile for component 'n', or
267     * null if all code-blocks for the current tile have been returned.
268     *
269     * @see CBlkWTData
270     * */
271    public final CBlkWTData getNextInternCodeBlock(int c, CBlkWTData cblk) {
272        // NOTE: this method is declared final since getNextCodeBlock() relies
273        // on this particular implementation
274        int k,j;
275        int tmp,shiftBits,jmin;
276        int w,h;
277        int outarr[];
278        float infarr[] = null;
279        CBlkWTDataFloat infblk;
280        float invstep;    // The inverse of the quantization step size
281        boolean intq;     // flag for quantizig ints
282        SubbandAn sb;
283        float stepUDR;    // The quantization step size (for a dynamic
284                          // range of 1, or unit)
285        int g = ((Integer)gbs.getTileCompVal(tIdx,c)).intValue();
286
287        // Are we quantizing ints or floats?
288        intq = (src.getDataType(tIdx,c) == DataBlk.TYPE_INT);
289
290        // Check that we have an output object
291        if (cblk == null) {
292            cblk = new CBlkWTDataInt();
293        }
294
295        // Cache input float code-block
296        infblk = this.infblk;
297
298        // Get data to quantize. When quantizing int data 'cblk' is used to
299        // get the data to quantize and to return the quantized data as well,
300        // that's why 'getNextCodeBlock()' is used. This can not be done when
301        // quantizing float data because of the different data types, that's
302        // why 'getNextInternCodeBlock()' is used in that case.
303        if (intq) { // Source data is int
304            cblk = src.getNextCodeBlock(c,cblk);
305            if (cblk == null) {
306                return null; // No more code-blocks in current tile for comp.
307            }
308            // Input and output arrays are the same (for "in place" quant.)
309            outarr = (int[])cblk.getData();
310        }
311        else { // Source data is float
312            // Can not use 'cblk' to get float data, use 'infblk'
313            infblk = (CBlkWTDataFloat) src.getNextInternCodeBlock(c,infblk);
314            if (infblk == null) {
315                // Release buffer from infblk: this enables to garbage collect
316                // the big buffer when we are done with last code-block of
317                // component.
318                this.infblk.setData(null);
319                return null; // No more code-blocks in current tile for comp.
320            }
321            this.infblk = infblk; // Save local cache
322            infarr = (float[])infblk.getData();
323            // Get output data array and check that there is memory to put the
324            // quantized coeffs in
325            outarr = (int[]) cblk.getData();
326            if (outarr == null || outarr.length < infblk.w*infblk.h) {
327                outarr = new int[infblk.w*infblk.h];
328                cblk.setData(outarr);
329            }
330            cblk.m = infblk.m;
331            cblk.n = infblk.n;
332            cblk.sb = infblk.sb;
333            cblk.ulx = infblk.ulx;
334            cblk.uly = infblk.uly;
335            cblk.w = infblk.w;
336            cblk.h = infblk.h;
337            cblk.wmseScaling = infblk.wmseScaling;
338            cblk.offset = 0;
339            cblk.scanw = cblk.w;
340        }
341
342        // Cache width, height and subband of code-block
343        w = cblk.w;
344        h = cblk.h;
345        sb = cblk.sb;
346
347        if(isReversible(tIdx,c)) { // Reversible only for int data
348            cblk.magbits = g-1+src.getNomRangeBits(c)+sb.anGainExp;
349            shiftBits = 31-cblk.magbits;
350
351            // Update the convertFactor field
352            cblk.convertFactor = (1<<shiftBits);
353
354            // Since we used getNextCodeBlock() to get the int data then
355            // 'offset' is 0 and 'scanw' is the width of the code-block The
356            // input and output arrays are the same (i.e. "in place")
357            for(j=w*h-1; j>=0; j--){
358                tmp = (outarr[j]<<shiftBits);
359                outarr[j] = ((tmp < 0) ? (1<<31)|(-tmp) : tmp);
360            }
361        }
362        else{ // Non-reversible, use step size
363            float baseStep =
364                ((Float)qsss.getTileCompVal(tIdx,c)).floatValue();
365
366            // Calculate magnitude bits and quantization step size
367            if(isDerived(tIdx,c)){
368                cblk.magbits = g-1+sb.level-
369                    (int)Math.floor(Math.log(baseStep)/log2);
370                stepUDR = baseStep/(1<<sb.level);
371            }
372            else{
373                cblk.magbits = g-1-(int)Math.floor(Math.log(baseStep/
374                                        (sb.l2Norm*(1<<sb.anGainExp)))/
375                                                   log2);
376                stepUDR = baseStep/(sb.l2Norm*(1<<sb.anGainExp));
377            }
378            shiftBits = 31-cblk.magbits;
379            // Calculate step that decoder will get and use that one.
380            stepUDR =
381                convertFromExpMantissa(convertToExpMantissa(stepUDR));
382            invstep = 1.0f/((1L<<(src.getNomRangeBits(c)+sb.anGainExp))*
383                            stepUDR);
384            // Normalize to magnitude bits (output fractional point)
385            invstep *= (1<<(shiftBits-src.getFixedPoint(c)));
386
387            // Update convertFactor and stepSize fields
388            cblk.convertFactor = invstep;
389            cblk.stepSize = ((1L<<(src.getNomRangeBits(c)+sb.anGainExp))*
390                             stepUDR);
391
392            if(intq){ // Quantizing int data
393                // Since we used getNextCodeBlock() to get the int data then
394                // 'offset' is 0 and 'scanw' is the width of the code-block
395                // The input and output arrays are the same (i.e. "in place")
396                for (j=w*h-1; j>=0; j--) {
397                    tmp = (int)(outarr[j]*invstep);
398                    outarr[j] = ((tmp < 0) ? (1<<31)|(-tmp) : tmp);
399                }
400            }
401            else { // Quantizing float data
402                for (j=w*h-1, k = infblk.offset+(h-1)*infblk.scanw+w-1,
403                         jmin = w*(h-1); j>=0; jmin -= w) {
404                    for (; j>=jmin; k--, j--) {
405                        tmp = (int)(infarr[k]*invstep);
406                        outarr[j] = ((tmp < 0) ? (1<<31)|(-tmp) : tmp);
407                    }
408                    // Jump to beggining of previous line in input
409                    k -= infblk.scanw - w;
410                }
411            }
412        }
413        // Return the quantized code-block
414        return cblk;
415    }
416
417    /**
418     * Calculates the parameters of the SubbandAn objects that depend on the
419     * Quantizer. The 'stepWMSE' field is calculated for each subband which is
420     * a leaf in the tree rooted at 'sb', for the specified component. The
421     * subband tree 'sb' must be the one for the component 'n'.
422     *
423     * @param sb The root of the subband tree.
424     *
425     * @param c The component index
426     *
427     * @see SubbandAn#stepWMSE
428     * */
429    protected void calcSbParams(SubbandAn sb,int c){
430        float baseStep;
431
432        if(sb.stepWMSE>0f) // parameters already calculated
433            return;
434        if(!sb.isNode){
435            if(isReversible(tIdx,c)){
436                sb.stepWMSE = (float) Math.pow(2,-(src.getNomRangeBits(c)<<1))*
437                    sb.l2Norm*sb.l2Norm;
438            }
439            else{
440                baseStep = ((Float)qsss.getTileCompVal(tIdx,c)).floatValue();
441                if(isDerived(tIdx,c)){
442                    sb.stepWMSE = baseStep*baseStep*
443                        (float)Math.pow(2,(sb.anGainExp-sb.level)<<1)*
444                        sb.l2Norm*sb.l2Norm;
445                }
446                else{
447                    sb.stepWMSE = baseStep*baseStep;
448                }
449            }
450        }
451        else{
452            calcSbParams((SubbandAn)sb.getLL(),c);
453            calcSbParams((SubbandAn)sb.getHL(),c);
454            calcSbParams((SubbandAn)sb.getLH(),c);
455            calcSbParams((SubbandAn)sb.getHH(),c);
456            sb.stepWMSE = 1f; // Signal that we already calculated this branch
457        }
458    }
459
460    /**
461     * Converts the floating point value to its exponent-mantissa
462     * representation. The mantissa occupies the 11 least significant bits
463     * (bits 10-0), and the exponent the previous 5 bits (bits 15-11).
464     *
465     * @param step The quantization step, normalized to a dynamic range of 1.
466     *
467     * @return The exponent mantissa representation of the step.
468     * */
469    public static int convertToExpMantissa(float step) {
470        int exp;
471
472        exp = (int)Math.ceil(-Math.log(step)/log2);
473        if (exp>QSTEP_MAX_EXPONENT) {
474            // If step size is too small for exponent representation, use the
475            // minimum, which is exponent QSTEP_MAX_EXPONENT and mantissa 0.
476            return (QSTEP_MAX_EXPONENT<<QSTEP_MANTISSA_BITS);
477        }
478        // NOTE: this formula does not support more than 5 bits for the
479        // exponent, otherwise (-1<<exp) might overflow (the - is used to be
480        // able to represent 2**31)
481        return (exp<<QSTEP_MANTISSA_BITS) |
482            ((int)((-step*(-1<<exp)-1f)*(1<<QSTEP_MANTISSA_BITS)+0.5f));
483    }
484
485    /**
486     * Converts the exponent-mantissa representation to its floating-point
487     * value. The mantissa occupies the 11 least significant bits (bits 10-0),
488     * and the exponent the previous 5 bits (bits 15-11).
489     *
490     * @param ems The exponent-mantissa representation of the step.
491     *
492     * @return The floating point representation of the step, normalized to a
493     * dynamic range of 1.
494     * */
495    private static float convertFromExpMantissa(int ems) {
496        // NOTE: this formula does not support more than 5 bits for the
497        // exponent, otherwise (-1<<exp) might overflow (the - is used to be
498        // able to represent 2**31)
499        return (-1f-((float)(ems&QSTEP_MAX_MANTISSA)) /
500                ((float)(1<<QSTEP_MANTISSA_BITS))) /
501            (float)(-1<<((ems>>QSTEP_MANTISSA_BITS)&QSTEP_MAX_EXPONENT));
502    }
503
504    /**
505     * Returns the maximum number of magnitude bits in any subband of the
506     * current tile.
507     *
508     * @param c the component number
509     *
510     * @return The maximum number of magnitude bits in all subbands of the
511     * current tile.
512     * */
513    public int getMaxMagBits(int c){
514        Subband sb = getAnSubbandTree(tIdx,c);
515        if(isReversible(tIdx,c)){
516            return getMaxMagBitsRev(sb,c);
517        }
518        else{
519            if(isDerived(tIdx,c)){
520                return getMaxMagBitsDerived(sb,tIdx,c);
521            }
522            else {
523                return getMaxMagBitsExpounded(sb,tIdx,c);
524            }
525        }
526    }
527
528
529    /**
530     * Returns the maximum number of magnitude bits in any subband of the
531     * current tile if reversible quantization is used
532     *
533     * @param sb The root of the subband tree of the current tile
534     *
535     * @param c the component number
536     *
537     * @return The highest number of magnitude bit-planes
538     * */
539    private int getMaxMagBitsRev(Subband sb, int c){
540        int tmp,max=0;
541        int g = ((Integer)gbs.getTileCompVal(tIdx,c)).intValue();
542
543        if(!sb.isNode)
544            return g-1+src.getNomRangeBits(c)+sb.anGainExp;
545
546        max=getMaxMagBitsRev(sb.getLL(),c);
547        tmp=getMaxMagBitsRev(sb.getLH(),c);
548        if(tmp>max)
549            max=tmp;
550        tmp=getMaxMagBitsRev(sb.getHL(),c);
551        if(tmp>max)
552            max=tmp;
553        tmp=getMaxMagBitsRev(sb.getHH(),c);
554        if(tmp>max)
555            max=tmp;
556
557        return max;
558    }
559
560    /**
561     * Returns the maximum number of magnitude bits in any subband in the
562     * given tile-component if derived quantization is used
563     *
564     * @param sb The root of the subband tree of the tile-component
565     *
566     * @param t Tile index
567     *
568     * @param c Component index
569     *
570     * @return The highest number of magnitude bit-planes
571     * */
572    private int getMaxMagBitsDerived(Subband sb,int t,int c){
573        int tmp,max=0;
574        int g = ((Integer)gbs.getTileCompVal(t,c)).intValue();
575
576        if(!sb.isNode){
577            float baseStep = ((Float)qsss.getTileCompVal(t,c)).floatValue();
578            return g-1+sb.level-(int)Math.floor(Math.log(baseStep)/log2);
579        }
580
581        max=getMaxMagBitsDerived(sb.getLL(),t,c);
582        tmp=getMaxMagBitsDerived(sb.getLH(),t,c);
583        if(tmp>max)
584            max=tmp;
585        tmp=getMaxMagBitsDerived(sb.getHL(),t,c);
586        if(tmp>max)
587            max=tmp;
588        tmp=getMaxMagBitsDerived(sb.getHH(),t,c);
589        if(tmp>max)
590            max=tmp;
591
592        return max;
593    }
594
595
596    /**
597     * Returns the maximum number of magnitude bits in any subband in the
598     * given tile-component if expounded quantization is used
599     *
600     * @param sb The root of the subband tree of the tile-component
601     *
602     * @param t Tile index
603     *
604     * @param c Component index
605     *
606     * @return The highest number of magnitude bit-planes
607     * */
608    private int getMaxMagBitsExpounded(Subband sb,int t,int c){
609        int tmp,max=0;
610        int g = ((Integer)gbs.getTileCompVal(t,c)).intValue();
611
612        if(!sb.isNode){
613            float baseStep = ((Float)qsss.getTileCompVal(t,c)).floatValue();
614            return g-1-
615                (int)Math.floor(Math.log(baseStep/
616                                (((SubbandAn)sb).l2Norm*(1<<sb.anGainExp)))/
617                                log2);
618        }
619
620        max=getMaxMagBitsExpounded(sb.getLL(),t,c);
621        tmp=getMaxMagBitsExpounded(sb.getLH(),t,c);
622        if(tmp>max)
623            max=tmp;
624        tmp=getMaxMagBitsExpounded(sb.getHL(),t,c);
625        if(tmp>max)
626            max=tmp;
627        tmp=getMaxMagBitsExpounded(sb.getHH(),t,c);
628        if(tmp>max)
629            max=tmp;
630
631        return max;
632    }
633}