Logo Search packages:      
Sourcecode: openexr version File versions  Download package

ImfWav.cpp

///////////////////////////////////////////////////////////////////////////
//
// Copyright (c) 2002, Industrial Light & Magic, a division of Lucas
// Digital Ltd. LLC
// 
// All rights reserved.
// 
// Redistribution and use in source and binary forms, with or without
// modification, are permitted provided that the following conditions are
// met:
// *       Redistributions of source code must retain the above copyright
// notice, this list of conditions and the following disclaimer.
// *       Redistributions in binary form must reproduce the above
// copyright notice, this list of conditions and the following disclaimer
// in the documentation and/or other materials provided with the
// distribution.
// *       Neither the name of Industrial Light & Magic nor the names of
// its contributors may be used to endorse or promote products derived
// from this software without specific prior written permission. 
// 
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
// "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
// LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
// A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
// OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
// SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
// LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
// THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
// (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
// OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
//
///////////////////////////////////////////////////////////////////////////



//-----------------------------------------------------------------------------
//
//    16-bit Haar Wavelet encoding and decoding
//
//    The source code in this file is derived from the encoding
//    and decoding routines written by Christian Rouet for his
//    PIZ image file format.
//
//-----------------------------------------------------------------------------


#include <ImfWav.h>

namespace Imf {
namespace {


//
// Wavelet basis functions without modulo arithmetic; they produce
// the best compression ratios when the wavelet-transformed data are
// Huffman-encoded, but the wavelet transform works only for 14-bit
// data (untransformed data values must be less than (1 << 14)).
//

inline void
wenc14 (unsigned short  a, unsigned short  b,
        unsigned short &l, unsigned short &h)
{
    short as = a;
    short bs = b;

    short ms = (as + bs) >> 1;
    short ds = as - bs;

    l = ms;
    h = ds;
}


inline void
wdec14 (unsigned short  l, unsigned short  h,
        unsigned short &a, unsigned short &b)
{
    short ls = l;
    short hs = h;

    int hi = hs;
    int ai = ls + (hi & 1) + (hi >> 1);

    short as = ai;
    short bs = ai - hi;

    a = as;
    b = bs;
}


//
// Wavelet basis functions with modulo arithmetic; they work with full
// 16-bit data, but Huffman-encoding the wavelet-transformed data doesn't
// compress the data quite as well.
//

const int NBITS = 16;
const int A_OFFSET =  1 << (NBITS  - 1);
const int M_OFFSET =  1 << (NBITS  - 1);
const int MOD_MASK = (1 <<  NBITS) - 1;


inline void
wenc16 (unsigned short  a, unsigned short  b,
        unsigned short &l, unsigned short &h)
{
    int ao =  (a + A_OFFSET) & MOD_MASK;
    int m  = ((ao + b) >> 1);
    int d  =   ao - b;

    if (d < 0)
      m = (m + M_OFFSET) & MOD_MASK;

    d &= MOD_MASK;

    l = m;
    h = d;
}


inline void
wdec16 (unsigned short  l, unsigned short  h,
        unsigned short &a, unsigned short &b)
{
    int m = l;
    int d = h;
    int bb = (m - (d >> 1)) & MOD_MASK;
    int aa = (d + bb - A_OFFSET) & MOD_MASK;
    b = bb;
    a = aa;
}

} // namespace


//
// 2D Wavelet encoding:
//

void
wav2Encode
    (unsigned short*    in,   // io: values are transformed in place
     int          nx,   // i : x size
     int          ox,   // i : x offset
     int          ny,   // i : y size
     int          oy,   // i : y offset
     unsigned short     mx)   // i : maximum in[x][y] value
{
    bool w14 = (mx < (1 << 14));
    int     n  = (nx > ny)? ny: nx;
    int     p  = 1;                 // == 1 <<  level
    int p2 = 2;               // == 1 << (level+1)

    //
    // Hierachical loop on smaller dimension n
    //

    while (p2 <= n)
    {
      unsigned short *py = in;
      unsigned short *ey = in + oy * (ny - p2);
      int oy1 = oy * p;
      int oy2 = oy * p2;
      int ox1 = ox * p;
      int ox2 = ox * p2;
      unsigned short i00,i01,i10,i11;

      //
      // Y loop
      //

      for (; py <= ey; py += oy2)
      {
          unsigned short *px = py;
          unsigned short *ex = py + ox * (nx - p2);

          //
          // X loop
          //

          for (; px <= ex; px += ox2)
          {
            unsigned short *p01 = px  + ox1;
            unsigned short *p10 = px  + oy1;
            unsigned short *p11 = p10 + ox1;

            //
            // 2D wavelet encoding
            //

            if (w14)
            {
                wenc14 (*px,  *p01, i00, i01);
                wenc14 (*p10, *p11, i10, i11);
                wenc14 (i00, i10, *px,  *p10);
                wenc14 (i01, i11, *p01, *p11);
            }
            else
            {
                wenc16 (*px,  *p01, i00, i01);
                wenc16 (*p10, *p11, i10, i11);
                wenc16 (i00, i10, *px,  *p10);
                wenc16 (i01, i11, *p01, *p11);
            }
          }

          //
          // Encode (1D) odd column (still in Y loop)
          //

          if (nx & p)
          {
            unsigned short *p10 = px + oy1;

            if (w14)
                wenc14 (*px, *p10, i00, *p10);
            else
                wenc16 (*px, *p10, i00, *p10);

            *px= i00;
          }
      }

      //
      // Encode (1D) odd line (must loop in X)
      //

      if (ny & p)
      {
          unsigned short *px = py;
          unsigned short *ex = py + ox * (nx - p2);

          for (; px <= ex; px += ox2)
          {
            unsigned short *p01 = px + ox1;

            if (w14)
                wenc14 (*px, *p01, i00, *p01);
            else
                wenc16 (*px, *p01, i00, *p01);

            *px= i00;
          }
      }

      //
      // Next level
      //

      p = p2;
      p2 <<= 1;
    }
}


//
// 2D Wavelet decoding:
//

void
wav2Decode
    (unsigned short*    in,   // io: values are transformed in place
     int          nx,   // i : x size
     int          ox,   // i : x offset
     int          ny,   // i : y size
     int          oy,   // i : y offset
     unsigned short     mx)   // i : maximum in[x][y] value
{
    bool w14 = (mx < (1 << 14));
    int     n = (nx > ny)? ny: nx;
    int     p = 1;
    int p2;

    //
    // Search max level
    //

    while (p <= n)
      p <<= 1;

    p >>= 1;
    p2 = p;
    p >>= 1;

    //
    // Hierarchical loop on smaller dimension n
    //

    while (p >= 1)
    {
      unsigned short *py = in;
      unsigned short *ey = in + oy * (ny - p2);
      int oy1 = oy * p;
      int oy2 = oy * p2;
      int ox1 = ox * p;
      int ox2 = ox * p2;
      unsigned short i00,i01,i10,i11;

      //
      // Y loop
      //

      for (; py <= ey; py += oy2)
      {
          unsigned short *px = py;
          unsigned short *ex = py + ox * (nx - p2);

          //
          // X loop
          //

          for (; px <= ex; px += ox2)
          {
            unsigned short *p01 = px  + ox1;
            unsigned short *p10 = px  + oy1;
            unsigned short *p11 = p10 + ox1;

            //
            // 2D wavelet decoding
            //

            if (w14)
            {
                wdec14 (*px,  *p10, i00, i10);
                wdec14 (*p01, *p11, i01, i11);
                wdec14 (i00, i01, *px,  *p01);
                wdec14 (i10, i11, *p10, *p11);
            }
            else
            {
                wdec16 (*px,  *p10, i00, i10);
                wdec16 (*p01, *p11, i01, i11);
                wdec16 (i00, i01, *px,  *p01);
                wdec16 (i10, i11, *p10, *p11);
            }
          }

          //
          // Decode (1D) odd column (still in Y loop)
          //

          if (nx & p)
          {
            unsigned short *p10 = px + oy1;

            if (w14)
                wdec14 (*px, *p10, i00, *p10);
            else
                wdec16 (*px, *p10, i00, *p10);

            *px= i00;
          }
      }

      //
      // Decode (1D) odd line (must loop in X)
      //

      if (ny & p)
      {
          unsigned short *px = py;
          unsigned short *ex = py + ox * (nx - p2);

          for (; px <= ex; px += ox2)
          {
            unsigned short *p01 = px + ox1;

            if (w14)
                wdec14 (*px, *p01, i00, *p01);
            else
                wdec16 (*px, *p01, i00, *p01);

            *px= i00;
          }
      }

      //
      // Next level
      //

      p2 = p;
      p >>= 1;
    }
}


} // namespace Imf

Generated by  Doxygen 1.6.0   Back to index