///////////////////////////////////////////////////////////////////////////

//

// Copyright (c) 2009-2014 DreamWorks Animation 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 DreamWorks Animation 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.

//

///////////////////////////////////////////////////////////////////////////

//---------------------------------------------------

//

// class DwaCompressor -- Store lossy RGB data by quantizing

//                          DCT components.

//

// First, we try and figure out what compression strategy to take

// based in channel name. For RGB channels, we want a lossy method

// described below. But, if we have alpha, we should do something

// different (and probably using RLE). If we have depth, or velocity,

// or something else, just fall back to ZIP. The rules for deciding 

// which strategy to use are setup in initializeDefaultChannelRules().

// When writing a file, the relevant rules needed to decode are written

// into the start of the data block, making a self-contained file. 

// If initializeDefaultChannelRules() doesn't quite suite your naming

// conventions, you can adjust the rules without breaking decoder

// compatability.

//

// If we're going to lossy compress R, G, or B channels, it's easier

// to toss bits in a more perceptual uniform space. One could argue

// at length as to what constitutes perceptually uniform, expecially 

// when storing either scene/input/focal plane referred and output referred

// data. 

//

// We'll compromise. For values <= 1, we use a traditional power function

// (without any of that straight-line business at the bottom). For values > 1,

// we want something more like a log function, since power functions blow

// up. At 1, we want a smooth blend between the functions. So, we use a 

// piecewise function that does just that - see dwaLookups.cpp for 

// a little more detail.

//

// Also, if we find that we have R, G, and B channels from the same layer,

// we can get a bit more compression efficiency by transforming to a Y'CbCr

// space. We use the 709 transform, but with Cb,Cr = 0 for an input of 

// (0, 0, 0), instead of the traditional Cb,Cr = .5. Shifting the zero point

// makes no sense with large range data. Transforms are done to from 

// the perceptual space data, not the linear-light space data (R'G'B' ->

// (Y'CbCr, not RGB -> YCbCr).

//

// Next, we forward DCT the data. This is done with a floating

// point DCT, as we don't really have control over the src range. The 

// resulting values are dropped to half-float precision. 

//

// Now, we need to quantize. Quantization departs from the usual way 

// of dividing and rounding. Instead, we start with some floating 

// point "base-error" value. From this, we can derive quantization 

// error for each DCT component. Take the standard JPEG quantization

// tables and normalize them by the smallest value. Then, multiply

// the normalized quant tables by our base-error value. This gives

// a range of errors for each DCT component.

//

// For each DCT component, we want to find a quantized value that 

// is within +- the per-component error. Pick the quantized value

// that has the fewest bits set in its' binary representation. 

// Brute-forcing the search would make for extremly inefficient 

// compression. Fortunatly, we can precompute a table to assist 

// with this search. 

//

// For each 16-bit float value, there are at most 15 other values with

// fewer bits set. We can precompute these values in a compact form, since

// many source values have far fewer that 15 possible quantized values. 

// Now, instead of searching the entire range +- the component error,

// we can just search at most 15 quantization candidates. The search can

// be accelerated a bit more by sorting the candidates by the 

// number of bits set, in increasing order. Then, the search can stop

// once a candidate is found w/i the per-component quantization 

// error range.

//

// The quantization strategy has the side-benefit that there is no

// de-quantization step upon decode, so we don't bother recording

// the quantization table.

//

// Ok. So we now have quantized values. Time for entropy coding. We

// can use either static Huffman or zlib/DEFLATE. The static Huffman

// is more efficient at compacting data, but can have a greater 

// overhead, especially for smaller tile/strip sizes. 

//

// There is some additional fun, like ZIP compressing the DC components

// instead of Huffman/zlib, which helps make things slightly smaller.

//

// Compression level is controlled by setting an int/float/double attribute

// on the header named "dwaCompressionLevel". This is a thinly veiled name for 

// the "base-error" value mentioned above. The "base-error" is just

// dwaCompressionLevel / 100000. The default value of 45.0 is generally 

// pretty good at generating "visually lossless" values at reasonable

// data rates. Setting dwaCompressionLevel to 0 should result in no additional

// quantization at the quantization stage (though there may be 

// quantization in practice at the CSC/DCT steps). But if you really

// want lossless compression, there are pleanty of other choices 

// of compressors ;)

//

// When dealing with FLOAT source buffers, we first quantize the source

// to HALF and continue down as we would for HALF source.

//

//---------------------------------------------------

#include "ImfDwaCompressor.h"

#include "ImfDwaCompressorSimd.h"

#include "ImfChannelList.h"

#include "ImfStandardAttributes.h"

#include "ImfHeader.h"

#include "ImfHuf.h"

#include "ImfInt64.h"

#include "ImfIntAttribute.h"

#include "ImfIO.h"

#include "ImfMisc.h"

#include "ImfNamespace.h"

#include "ImfRle.h"

#include "ImfSimd.h"

#include "ImfSystemSpecific.h"

#include "ImfXdr.h"

#include "ImfZip.h"

#include "ImathFun.h"

#include "ImathBox.h"

#include "ImathVec.h"

#include "half.h"

#include "halfLimits.h"

#include "dwaLookups.h"

#include <vector>

#include <string>

#include <cctype>

#include <cassert>

#include <algorithm>

// Windows specific addition to prevent the indirect import of the redefined min/max macros

#if defined _WIN32 || defined _WIN64

  #ifdef NOMINMAX

    #undef NOMINMAX

  #endif

  #define NOMINMAX

#endif

#include <zlib.h>

OPENEXR_IMF_INTERNAL_NAMESPACE_SOURCE_ENTER

namespace {

    //

    // Function pointer to dispatch to an approprate 

    // convertFloatToHalf64_* impl, based on runtime cpu checking.

    // Should be initialized in DwaCompressor::initializeFuncs()

    //

    void (*convertFloatToHalf64)(unsigned short*, float*) =

        convertFloatToHalf64_scalar;

    // 

    // Function pointer for dispatching a fromHalfZigZag_ impl

    //

    void (*fromHalfZigZag)(unsigned short*, float*) =

        fromHalfZigZag_scalar;

    //

    // Dispatch the inverse DCT on an 8x8 block, where the last

    // n rows can be all zeros. The n=0 case converts the full block.

    //

    void (*dctInverse8x8_0)(float*) = dctInverse8x8_scalar<0>;

    void (*dctInverse8x8_1)(float*) = dctInverse8x8_scalar<1>;

    void (*dctInverse8x8_2)(float*) = dctInverse8x8_scalar<2>;

    void (*dctInverse8x8_3)(float*) = dctInverse8x8_scalar<3>;

    void (*dctInverse8x8_4)(float*) = dctInverse8x8_scalar<4>;

    void (*dctInverse8x8_5)(float*) = dctInverse8x8_scalar<5>;

    void (*dctInverse8x8_6)(float*) = dctInverse8x8_scalar<6>;

    void (*dctInverse8x8_7)(float*) = dctInverse8x8_scalar<7>;

} // namespace

struct DwaCompressor::ChannelData

{

    std::string         name;

    CompressorScheme    compression;  

    int                 xSampling;

    int                 ySampling;

    PixelType           type;

    bool                pLinear;

    int                 width;

    int                 height;

    //

    // Incoming and outgoing data is scanline interleaved, and it's much

    // easier to operate on contiguous data.  Assuming the planare unc

    // buffer is to hold RLE data, we need to rearrange to make bytes

    // adjacent.

    //

    char               *planarUncBuffer;

    char               *planarUncBufferEnd;

    char               *planarUncRle[4];

    char               *planarUncRleEnd[4];

    PixelType           planarUncType;

    int                 planarUncSize;

};

struct DwaCompressor::CscChannelSet

{

    int idx[3];

};

struct DwaCompressor::Classifier

{

    Classifier (std::string suffix,

                CompressorScheme scheme,

                PixelType type,

                int cscIdx,

                bool caseInsensitive):

        _suffix(suffix),

        _scheme(scheme),

        _type(type),

        _cscIdx(cscIdx),

        _caseInsensitive(caseInsensitive)

    {

        if (caseInsensitive) 

            std::transform(_suffix.begin(), _suffix.end(), _suffix.begin(), tolower);

    }

    Classifier (const char *&ptr, int size)

    {

        if (size <= 0) 

            throw IEX_NAMESPACE::InputExc("Error uncompressing DWA data"

                                " (truncated rule).");

        {

            char suffix[Name::SIZE];

            memset (suffix, 0, Name::SIZE);

            Xdr::read<CharPtrIO> (ptr, std::min(size, Name::SIZE-1), suffix);

            _suffix = std::string(suffix);

        }

        if (size < _suffix.length() + 1 + 2*Xdr::size<char>()) 

            throw IEX_NAMESPACE::InputExc("Error uncompressing DWA data"

                                " (truncated rule).");

        char value;

        Xdr::read<CharPtrIO> (ptr, value);

        _cscIdx = (int)(value >> 4) - 1;

        if (_cscIdx < -1 || _cscIdx >= 3) 

            throw IEX_NAMESPACE::InputExc("Error uncompressing DWA data"

                                " (corrupt cscIdx rule).");

        _scheme = (CompressorScheme)((value >> 2) & 3);

        if (_scheme < 0 || _scheme >= NUM_COMPRESSOR_SCHEMES) 

            throw IEX_NAMESPACE::InputExc("Error uncompressing DWA data"

                                " (corrupt scheme rule).");

        _caseInsensitive = (value & 1 ? true : false);

        Xdr::read<CharPtrIO> (ptr, value);

        if (value < 0 || value >= NUM_PIXELTYPES) 

            throw IEX_NAMESPACE::InputExc("Error uncompressing DWA data"

                                " (corrupt rule).");

        _type = (PixelType)value;

    }

    bool match (const std::string &suffix, const PixelType type) const

    {

        if (_type != type) return false;

        if (_caseInsensitive) 

        {

            std::string tmp(suffix);

            std::transform(tmp.begin(), tmp.end(), tmp.begin(), tolower);

            return tmp == _suffix;

        }

        return suffix == _suffix;

    }

    size_t size () const 

    {

        // string length + \0

        size_t sizeBytes = _suffix.length() + 1;

        // 1 byte for scheme / cscIdx / caseInsensitive, and 1 byte for type

        sizeBytes += 2 * Xdr::size<char>();

        return sizeBytes;

    }

    void write (char *&ptr) const

    {

        Xdr::write<CharPtrIO> (ptr, _suffix.c_str());

        // Encode _cscIdx (-1-3) in the upper 4 bits,

        //        _scheme (0-2)  in the next 2 bits

        //        _caseInsen     in the bottom bit

        unsigned char value = 0;

        value |= ((unsigned char)(_cscIdx+1)      & 15) << 4;

        value |= ((unsigned char)_scheme          &  3) << 2;

        value |=  (unsigned char)_caseInsensitive &  1;

        Xdr::write<CharPtrIO> (ptr, value);

        Xdr::write<CharPtrIO> (ptr, (unsigned char)_type);

    }

    std::string      _suffix;

    CompressorScheme _scheme;

    PixelType        _type;

    int              _cscIdx;

    bool             _caseInsensitive;

};

//

// Base class for the LOSSY_DCT decoder classes

//

class DwaCompressor::LossyDctDecoderBase

{

  public:

    LossyDctDecoderBase

        (char *packedAc,

         char *packedDc,

         const unsigned short *toLinear,

         int width,

         int height);

    virtual ~LossyDctDecoderBase ();

    void execute();

    //

    // These return number of items, not bytes. Each item

    // is an unsigned short

    //

    int numAcValuesEncoded() const { return _packedAcCount; }

    int numDcValuesEncoded() const { return _packedDcCount; }

  protected:

    //

    // Un-RLE the packed AC components into 

    // a half buffer. The half block should 

    // be the full 8x8 block (in zig-zag order

    // still), not the first AC component. 

    //

    // currAcComp is advanced as bytes are decoded.

    //

    // This returns the index of the last non-zero

    // value in the buffer - with the index into zig zag

    // order data. If we return 0, we have DC only data.

    // 

    int unRleAc (unsigned short *&currAcComp,

                 unsigned short  *halfZigBlock); 

    //

    // if NATIVE and XDR are really the same values, we can

    // skip some processing and speed things along

    //

    bool                  _isNativeXdr;

    //

    // Counts of how many items have been packed into the

    // AC and DC buffers

    //

    int                   _packedAcCount;

    int                   _packedDcCount;

    //

    // AC and DC buffers to pack

    //

    char                 *_packedAc;

    char                 *_packedDc;

    // 

    // half -> half LUT to transform from nonlinear to linear

    //

    const unsigned short *_toLinear;

    //

    // image dimensions

    //

    int                   _width;

    int                   _height;

    //

    // Pointers to the start of each scanlines, to be filled on decode

    // Generally, these will be filled by the subclasses.

    //

    std::vector< std::vector<char *> > _rowPtrs;

    // 

    // The type of each data that _rowPtrs[i] is referring. Layout

    // is in the same order as _rowPtrs[].

    //

    std::vector<PixelType>             _type;

    std::vector<SimdAlignedBuffer64f>  _dctData;

};

//

// Used to decode a single channel of LOSSY_DCT data.

//

class DwaCompressor::LossyDctDecoder: public LossyDctDecoderBase

{

  public:

    //

    // toLinear is a half-float LUT to convert the encoded values 

    // back to linear light. If you want to skip this step, pass

    // in NULL here.

    //

    LossyDctDecoder

        (std::vector<char *> &rowPtrs,

         char *packedAc,

         char *packedDc,

         const unsigned short *toLinear,

         int width,

         int height,

         PixelType type)

    :

        LossyDctDecoderBase(packedAc, packedDc, toLinear, width, height)

    {

        _rowPtrs.push_back(rowPtrs);

        _type.push_back(type);

    }

    virtual ~LossyDctDecoder () {}

};

//

// Used to decode 3 channels of LOSSY_DCT data that

// are grouped together and color space converted.

//

class DwaCompressor::LossyDctDecoderCsc: public LossyDctDecoderBase

{

  public:

    //

    // toLinear is a half-float LUT to convert the encoded values 

    // back to linear light. If you want to skip this step, pass

    // in NULL here.

    //

    LossyDctDecoderCsc

        (std::vector<char *> &rowPtrsR,

         std::vector<char *> &rowPtrsG,

         std::vector<char *> &rowPtrsB,

         char *packedAc,

         char *packedDc,

         const unsigned short *toLinear,

         int width,

         int height,

         PixelType typeR,

         PixelType typeG,

         PixelType typeB)

    :

        LossyDctDecoderBase(packedAc, packedDc, toLinear, width, height)

    {

        _rowPtrs.push_back(rowPtrsR);

        _rowPtrs.push_back(rowPtrsG);

        _rowPtrs.push_back(rowPtrsB);

        _type.push_back(typeR);

        _type.push_back(typeG);

        _type.push_back(typeB);

    }

    virtual ~LossyDctDecoderCsc () {}

};

// 

// Base class for encoding using the lossy DCT scheme

//

class DwaCompressor::LossyDctEncoderBase

{

  public:

    LossyDctEncoderBase

        (float quantBaseError,

         char *packedAc,

         char *packedDc,

         const unsigned short *toNonlinear,

         int width,

         int height);

    virtual ~LossyDctEncoderBase ();

    void execute ();

    //

    // These return number of items, not bytes. Each item

    // is an unsigned short

    //

    int     numAcValuesEncoded () const {return _numAcComp;}

    int     numDcValuesEncoded () const {return _numDcComp;}

  protected:

    void    toZigZag (half *dst, half *src);

    int     countSetBits (unsigned short src);

    half    quantize (half src, float errorTolerance);

    void    rleAc (half *block, unsigned short *&acPtr);

    float                      _quantBaseError;

    int                        _width,

                               _height;

    const unsigned short      *_toNonlinear;

    int                        _numAcComp,

                               _numDcComp;

    std::vector< std::vector<const char *> > _rowPtrs;

    std::vector<PixelType>                   _type;

    std::vector<SimdAlignedBuffer64f>        _dctData;

    //

    // Pointers to the buffers where AC and DC

    // DCT components should be packed for 

    // lossless compression downstream

    //

    char                      *_packedAc;

    char                      *_packedDc;

    //

    // Our "quantization tables" - the example JPEG tables, 

    // normalized so that the smallest value in each is 1.0.

    // This gives us a relationship between error in DCT 

    // components

    //

    float                      _quantTableY[64];

    float                      _quantTableCbCr[64];

};

//

// Single channel lossy DCT encoder

//

class DwaCompressor::LossyDctEncoder: public LossyDctEncoderBase

{

  public:

    LossyDctEncoder

        (float quantBaseError,

         std::vector<const char *> &rowPtrs,

         char *packedAc,

         char *packedDc,

         const unsigned short *toNonlinear,

         int width,

         int height,

         PixelType type)

    :

        LossyDctEncoderBase

            (quantBaseError, packedAc, packedDc, toNonlinear, width, height)

    {

        _rowPtrs.push_back(rowPtrs);

        _type.push_back(type);

    }

    virtual ~LossyDctEncoder () {}

};

//

// RGB channel lossy DCT encoder

//

class DwaCompressor::LossyDctEncoderCsc: public LossyDctEncoderBase

{

  public:

    LossyDctEncoderCsc

        (float quantBaseError,

         std::vector<const char *> &rowPtrsR,

         std::vector<const char *> &rowPtrsG,

         std::vector<const char *> &rowPtrsB,

         char *packedAc,

         char *packedDc,

         const unsigned short *toNonlinear,

         int width,

         int height,

         PixelType typeR,

         PixelType typeG,

         PixelType typeB)

    :

        LossyDctEncoderBase

            (quantBaseError, packedAc, packedDc, toNonlinear, width, height)

    {

        _type.push_back(typeR);

        _type.push_back(typeG);

        _type.push_back(typeB);

        _rowPtrs.push_back(rowPtrsR);

        _rowPtrs.push_back(rowPtrsG);

        _rowPtrs.push_back(rowPtrsB);

    }

    virtual ~LossyDctEncoderCsc () {}

};

// ==============================================================

//

//                     LossyDctDecoderBase

//

// --------------------------------------------------------------

DwaCompressor::LossyDctDecoderBase::LossyDctDecoderBase

    (char *packedAc,

     char *packedDc,

     const unsigned short *toLinear,

     int width,

     int height)

:

    _isNativeXdr(false),

    _packedAcCount(0),

    _packedDcCount(0),

    _packedAc(packedAc),

    _packedDc(packedDc),

    _toLinear(toLinear),

    _width(width),

    _height(height)

{

    if (_toLinear == 0)

        _toLinear = get_dwaCompressorNoOp();

    _isNativeXdr = GLOBAL_SYSTEM_LITTLE_ENDIAN;

}

DwaCompressor::LossyDctDecoderBase::~LossyDctDecoderBase () {}

void

DwaCompressor::LossyDctDecoderBase::execute ()

{

    int numComp        = _rowPtrs.size();

    int lastNonZero    = 0;

    int numBlocksX     = (int) ceil ((float)_width  / 8.0f);

    int numBlocksY     = (int) ceil ((float)_height / 8.0f);

    int leftoverX      = _width  - (numBlocksX-1) * 8;

    int leftoverY      = _height - (numBlocksY-1) * 8;

    int numFullBlocksX = (int)floor ((float)_width / 8.0f);

    unsigned short tmpShortNative = 0;

    unsigned short tmpShortXdr    = 0;

    const char *tmpConstCharPtr   = 0;

    unsigned short                    *currAcComp = (unsigned short *)_packedAc;

    std::vector<unsigned short *>      currDcComp (_rowPtrs.size());

    std::vector<SimdAlignedBuffer64us> halfZigBlock (_rowPtrs.size());

    if (_type.size() != _rowPtrs.size())

        throw IEX_NAMESPACE::BaseExc ("Row pointers and types mismatch in count");

    if ((_rowPtrs.size() != 3) && (_rowPtrs.size() != 1))

        throw IEX_NAMESPACE::NoImplExc ("Only 1 and 3 channel encoding is supported");

    _dctData.resize(numComp);

    //

    // Allocate a temp aligned buffer to hold a rows worth of full 

    // 8x8 half-float blocks

    //

    unsigned char *rowBlockHandle = new unsigned char

        [numComp * numBlocksX * 64 * sizeof(unsigned short) + _SSE_ALIGNMENT];

    unsigned short *rowBlock[3];

    rowBlock[0] = (unsigned short*)rowBlockHandle;

    for (int i = 0; i < _SSE_ALIGNMENT; ++i)

    {

        if (((size_t)(rowBlockHandle + i) & _SSE_ALIGNMENT_MASK) == 0)

            rowBlock[0] = (unsigned short *)(rowBlockHandle + i);

    }

    for (int comp = 1; comp < numComp; ++comp)

        rowBlock[comp] = rowBlock[comp - 1] + numBlocksX * 64;

    //

    // Pack DC components together by common plane, so we can get 

    // a little more out of differencing them. We'll always have 

    // one component per block, so we can computed offsets.

    //

    currDcComp[0] = (unsigned short *)_packedDc;

    for (unsigned int comp = 1; comp < numComp; ++comp)

        currDcComp[comp] = currDcComp[comp - 1] + numBlocksX * numBlocksY;

    for (int blocky = 0; blocky < numBlocksY; ++blocky)

    {

        int maxY = 8;

        if (blocky == numBlocksY-1)

            maxY = leftoverY;

        int maxX = 8;

        for (int blockx = 0; blockx < numBlocksX; ++blockx)

        {

            if (blockx == numBlocksX-1)

                maxX = leftoverX;

            //

            // If we can detect that the block is constant values

            // (all components only have DC values, and all AC is 0),

            // we can do everything only on 1 value, instead of all

            // 64. 

            //

            // This won't really help for regular images, but it is

            // meant more for layers with large swaths of black 

            //

            bool blockIsConstant = true;

            for (unsigned int comp = 0; comp < numComp; ++comp)

            {

                //

                // DC component is stored separately

                //

                #ifdef IMF_HAVE_SSE2

                    {

                        __m128i *dst = (__m128i*)halfZigBlock[comp]._buffer;

                        dst[7] = _mm_setzero_si128();

                        dst[6] = _mm_setzero_si128();

                        dst[5] = _mm_setzero_si128();

                        dst[4] = _mm_setzero_si128();

                        dst[3] = _mm_setzero_si128();

                        dst[2] = _mm_setzero_si128();

                        dst[1] = _mm_setzero_si128();

                        dst[0] = _mm_insert_epi16

                            (_mm_setzero_si128(), *currDcComp[comp]++, 0);

                    }

                #else  /* IMF_HAVE_SSE2 */

                    memset (halfZigBlock[comp]._buffer, 0, 64 * 2);

                    halfZigBlock[comp]._buffer[0] = *currDcComp[comp]++;

                #endif /* IMF_HAVE_SSE2 */

                _packedDcCount++;

                //

                // UnRLE the AC. This will modify currAcComp

                //

                lastNonZero = unRleAc (currAcComp, halfZigBlock[comp]._buffer);

                //

                // Convert from XDR to NATIVE

                //

                if (!_isNativeXdr)

                {

                    for (int i = 0; i < 64; ++i)

                    {

                        tmpShortXdr      = halfZigBlock[comp]._buffer[i];

                        tmpConstCharPtr  = (const char *)&tmpShortXdr;

                        Xdr::read<CharPtrIO> (tmpConstCharPtr, tmpShortNative);

                        halfZigBlock[comp]._buffer[i] = tmpShortNative;

                    }

                }

                if (lastNonZero == 0)

                {

                    //

                    // DC only case - AC components are all 0   

                    //

                    half h;

                    h.setBits (halfZigBlock[comp]._buffer[0]);

                    _dctData[comp]._buffer[0] = (float)h;

                    dctInverse8x8DcOnly (_dctData[comp]._buffer);

                }

                else

                {

                    //

                    // We have some AC components that are non-zero. 

                    // Can't use the 'constant block' optimization

                    //

                    blockIsConstant = false;

                    //

                    // Un-Zig zag 

                    //

                    (*fromHalfZigZag)

                        (halfZigBlock[comp]._buffer, _dctData[comp]._buffer);

                    //

                    // Zig-Zag indices in normal layout are as follows:

                    //

                    // 0   1   5   6   14  15  27  28

                    // 2   4   7   13  16  26  29  42

                    // 3   8   12  17  25  30  41  43

                    // 9   11  18  24  31  40  44  53

                    // 10  19  23  32  39  45  52  54

                    // 20  22  33  38  46  51  55  60

                    // 21  34  37  47  50  56  59  61

                    // 35  36  48  49  57  58  62  63

                    //

                    // If lastNonZero is less than the first item on

                    // each row, we know that the whole row is zero and 

                    // can be skipped in the row-oriented part of the

                    // iDCT. 

                    //

                    // The unrolled logic here is:

                    //

                    //    if lastNonZero < rowStartIdx[i],

                    //    zeroedRows = rowsEmpty[i]

                    //

                    // where:

                    //

                    //    const int rowStartIdx[] = {2, 3, 9, 10, 20, 21, 35};

                    //    const int rowsEmpty[]   = {7, 6, 5,  4,  3,  2,  1};

                    //

                    if (lastNonZero < 2)

                        dctInverse8x8_7(_dctData[comp]._buffer);

                    else if (lastNonZero < 3)

                        dctInverse8x8_6(_dctData[comp]._buffer);

                    else if (lastNonZero < 9)

                        dctInverse8x8_5(_dctData[comp]._buffer);

                    else if (lastNonZero < 10)

                        dctInverse8x8_4(_dctData[comp]._buffer);

                    else if (lastNonZero < 20)

                        dctInverse8x8_3(_dctData[comp]._buffer);

                    else if (lastNonZero < 21)

                        dctInverse8x8_2(_dctData[comp]._buffer);

                    else if (lastNonZero < 35)

                        dctInverse8x8_1(_dctData[comp]._buffer);

                    else

                        dctInverse8x8_0(_dctData[comp]._buffer);

                }

            }

            //

            // Perform the CSC

            //

            if (numComp == 3)

            {

                if (!blockIsConstant)

                {

                    csc709Inverse64 (_dctData[0]._buffer, 

                                     _dctData[1]._buffer, 

                                     _dctData[2]._buffer);

                }

                else

                {

                    csc709Inverse (_dctData[0]._buffer[0], 

                                   _dctData[1]._buffer[0], 

                                   _dctData[2]._buffer[0]);

                }

            }

            //

            // Float -> Half conversion. 

            //

            // If the block has a constant value, just convert the first pixel.

            //

            for (unsigned int comp = 0; comp < numComp; ++comp)

            {

                if (!blockIsConstant)

                {

                    (*convertFloatToHalf64)

                        (&rowBlock[comp][blockx*64], _dctData[comp]._buffer);

                }

                else

                {

                    #ifdef IMF_HAVE_SSE2

                        __m128i *dst = (__m128i*)&rowBlock[comp][blockx*64];

                        dst[0] = _mm_set1_epi16

                            (((half)_dctData[comp]._buffer[0]).bits());

                        dst[1] = dst[0];

                        dst[2] = dst[0];

                        dst[3] = dst[0];

                        dst[4] = dst[0];

                        dst[5] = dst[0];

                        dst[6] = dst[0];

                        dst[7] = dst[0];

                    #else  /* IMF_HAVE_SSE2 */

                        unsigned short *dst = &rowBlock[comp][blockx*64];

                        dst[0] = ((half)_dctData[comp]._buffer[0]).bits();

                        for (int i = 1; i < 64; ++i)

                        {

                            dst[i] = dst[0];

                        }

                    #endif /* IMF_HAVE_SSE2 */

                } // blockIsConstant

            } // comp

        } // blockx

        //

        // At this point, we have half-float nonlinear value blocked

        // in rowBlock[][]. We need to unblock the data, transfer

        // back to linear, and write the results in the _rowPtrs[].

        //

        // There is a fast-path for aligned rows, which helps

        // things a little. Since this fast path is only valid

        // for full 8-element wide blocks, the partial x blocks

        // are broken into a separate loop below.

        //

        // At the moment, the fast path requires:

        //   * sse support

        //   * aligned row pointers

        //   * full 8-element wide blocks

        //

        for (int comp = 0; comp < numComp; ++comp)

        {

            //

            // Test if we can use the fast path

            //

        #ifdef IMF_HAVE_SSE2

            bool fastPath = true;

            for (int y = 8 * blocky; y < 8 * blocky + maxY; ++y)

            {

                if ((size_t)_rowPtrs[comp][y] & _SSE_ALIGNMENT_MASK)

                    fastPath = false;

            }

            if (fastPath)

            {

                //

                // Handle all the full X blocks, in a fast path with sse2 and

                // aligned row pointers

                //

                for (int y=8*blocky; y<8*blocky+maxY; ++y)

                {

                    __m128i *dst = (__m128i *)_rowPtrs[comp][y];

                    __m128i *src = (__m128i *)&rowBlock[comp][(y & 0x7) * 8];