/************************************************************************/
/*                          UncompressBlock()                           */
/*                                                                      */
/*      Uncompress ESRI Grid compression format block.                  */
/************************************************************************/

static CPLErr UncompressBlock( GByte *pabyCData, int /* nSrcBytes */,
                               GByte *pabyDest, int nMaxPixels, 
                               int nDataType )

{
    GUInt32  nDataMin, nDataOffset;
    int      nNumBits, nPixelsOutput=0;			
    GInt32   nNumRuns;
    GByte *pabyCounter, *pabyValues;
    int   nValueBitOffset;

    memcpy( &nDataMin, pabyCData, 4 );
    nDataMin = CPL_LSBWORD32( nDataMin );
        
    memcpy( &nNumRuns, pabyCData+4, 4 );
    nNumRuns = CPL_LSBWORD32( nNumRuns );
        
    memcpy( &nDataOffset, pabyCData+8, 4 );
    nDataOffset = CPL_LSBWORD32( nDataOffset );

    nNumBits = pabyCData[12];

/* ==================================================================== */
/*      If this is not run length encoded, but just reduced             */
/*      precision, handle it now.                                       */
/* ==================================================================== */
    if( nNumRuns == -1 )
    {
        pabyValues = pabyCData + 13;
        nValueBitOffset = 0;

/* -------------------------------------------------------------------- */
/*      For floating point output, the saved 16-bit value will be       */
/*      based on the minimum value for the tile. Note that the there    */
/*      still seems to be some problems selecting the correct factor    */
/*      right near a power of 2 boundary (like 512). I have not been    */
/*      able to determine how the switchover occurs at power of 2       */
/*      boundaries for the minimum value.                               */
/* -------------------------------------------------------------------- */
        float fMultFactor = 1.0;
        float fMinValue = *((float *) &nDataMin);
        if ( nDataType == EPT_f32 )
        {
            int nDataMin = ( fMinValue > 0.0 ) ? (int)fMinValue : (int)-fMinValue;
            int nDivShift = -9;
            for ( ; ( nDataMin > 0 ); nDivShift++, nDataMin >>= 1 ) {}
            if ( nDivShift < 0 )
            {
                nDivShift = -nDivShift;
                fMultFactor = 1.0 / ( 1 << nDivShift );
            }
            else
            {
                fMultFactor = ( 1 << nDivShift );
            }
        }

        for( nPixelsOutput = 0; nPixelsOutput < nMaxPixels; nPixelsOutput++ )
        {
            int	nDataValue, nRawValue;

/* -------------------------------------------------------------------- */
/*      Extract the data value in a way that depends on the number      */
/*      of bits in it.                                                  */
/* -------------------------------------------------------------------- */
            if( nNumBits == 0 )
            {
                nRawValue = 0;
            }
            else if( nNumBits == 1 )
            {
                nRawValue =
                    (pabyValues[nValueBitOffset>>3] >> (nValueBitOffset&7)) & 0x1;
                nValueBitOffset++;
            }
            else if( nNumBits == 2 )
            {
                nRawValue =
                    (pabyValues[nValueBitOffset>>3] >> (nValueBitOffset&7)) & 0x3;
                nValueBitOffset += 2;
            }
            else if( nNumBits == 4 )
            {
                nRawValue =
                    (pabyValues[nValueBitOffset>>3] >> (nValueBitOffset&7)) & 0xf;
                nValueBitOffset += 4;
            }
            else if( nNumBits == 8 )
            {
                nRawValue = *pabyValues;
                pabyValues++;
            }
            else if( nNumBits == 16 )
            {
                nRawValue = 256 * *(pabyValues++);
                nRawValue += *(pabyValues++);
            }
            else if( nNumBits == 32 )
            {
                nRawValue = 256 * 256 * 256 * *(pabyValues++);
                nRawValue += 256 * 256 * *(pabyValues++);
                nRawValue += 256 * *(pabyValues++);
                nRawValue += *(pabyValues++);
            }
            else
            {
                printf( "nNumBits = %d\n", nNumBits );
                CPLAssert( FALSE );
                nRawValue = 0;
            }

/* -------------------------------------------------------------------- */
/*      Offset by the minimum value.                                    */
/* -------------------------------------------------------------------- */
            nDataValue = nRawValue + nDataMin;

/* -------------------------------------------------------------------- */
/*      Now apply to the output buffer in a type specific way.          */
/* -------------------------------------------------------------------- */
            if( nDataType == EPT_u8 )
            {
                ((GByte *) pabyDest)[nPixelsOutput] = (GByte) nDataValue;
            }
            else if( nDataType == EPT_u1 )
            {
                if( nDataValue == 1 )
                    pabyDest[nPixelsOutput>>3] |= (1 << (nPixelsOutput & 0x7));
                else
                    pabyDest[nPixelsOutput>>3] &= ~(1<<(nPixelsOutput & 0x7));
            }
            else if( nDataType == EPT_u2 )
            {
                if( (nPixelsOutput & 0x1) == 0 )
                    pabyDest[nPixelsOutput>>2] = (GByte) nDataValue;
                else if( (nPixelsOutput & 0x1) == 1 )
                    pabyDest[nPixelsOutput>>2] |= (GByte) (nDataValue<<2);
                else if( (nPixelsOutput & 0x1) == 2 )
                    pabyDest[nPixelsOutput>>2] |= (GByte) (nDataValue<<4);
                else
                    pabyDest[nPixelsOutput>>2] |= (GByte) (nDataValue<<6);
            }
            else if( nDataType == EPT_u4 )
            {
                if( (nPixelsOutput & 0x1) == 0 )
                    pabyDest[nPixelsOutput>>1] = (GByte) nDataValue;
                else
                    pabyDest[nPixelsOutput>>1] |= (GByte) (nDataValue<<4);
            }
            else if( nDataType == EPT_u16 )
            {
                ((GUInt16 *) pabyDest)[nPixelsOutput] = (GUInt16) nDataValue;
            }
            else if( nDataType == EPT_s16 )
            {
                ((GInt16 *) pabyDest)[nPixelsOutput] = (GInt16) nDataValue;
            }
            else if( nDataType == EPT_s32 )
            {
                ((GInt32 *) pabyDest)[nPixelsOutput] = nDataValue;
            }
            else if( nDataType == EPT_u32 )
            {
                ((GUInt32 *) pabyDest)[nPixelsOutput] = nDataValue;
            }
/* -------------------------------------------------------------------- */
/*      Note, floating point values are handled somewhat                */
/*      differently, and I've only been able to test f32 with a         */
/*      16bit offset value (see bug #1000 and                           */
/*      data/imagine/bug1000/float.img)                                 */
/* -------------------------------------------------------------------- */
            else if( nDataType == EPT_f32 )
            {
                float fValue;

                if( nNumBits == 16 )
                    fValue = fMinValue + fMultFactor * (nRawValue / 32768.0);
                else
                    CPLAssert( FALSE );
                
                ((float *) pabyDest)[nPixelsOutput] = fValue;
            }
            else
            {
                CPLAssert( FALSE );
            }
        }

        return CE_None;
    }

/* ==================================================================== */
/*      Establish data pointers for runs.                               */
/* ==================================================================== */
    pabyCounter = pabyCData + 13;
    pabyValues = pabyCData + nDataOffset;
    nValueBitOffset = 0;
    
/* -------------------------------------------------------------------- */
/*      Loop over runs.                                                 */
/* -------------------------------------------------------------------- */
    int    iRun;

    for( iRun = 0; iRun < nNumRuns; iRun++ )
    {
        int	nRepeatCount = 0;
        int	nDataValue;

/* -------------------------------------------------------------------- */
/*      Get the repeat count.  This can be stored as one, two, three    */
/*      or four bytes depending on the low order two bits of the        */
/*      first byte.                                                     */
/* -------------------------------------------------------------------- */
        if( ((*pabyCounter) & 0xc0) == 0x00 )
        {
            nRepeatCount = (*(pabyCounter++)) & 0x3f;
        }
        else if( ((*pabyCounter) & 0xc0) == 0x40 )
        {
            nRepeatCount = (*(pabyCounter++)) & 0x3f;
            nRepeatCount = nRepeatCount * 256 + (*(pabyCounter++));
        }
        else if( ((*pabyCounter) & 0xc0) == 0x80 )
        {
            nRepeatCount = (*(pabyCounter++)) & 0x3f;
            nRepeatCount = nRepeatCount * 256 + (*(pabyCounter++));
            nRepeatCount = nRepeatCount * 256 + (*(pabyCounter++));
        }
        else if( ((*pabyCounter) & 0xc0) == 0xc0 )
        {
            nRepeatCount = (*(pabyCounter++)) & 0x3f;
            nRepeatCount = nRepeatCount * 256 + (*(pabyCounter++));
            nRepeatCount = nRepeatCount * 256 + (*(pabyCounter++));
            nRepeatCount = nRepeatCount * 256 + (*(pabyCounter++));
        }

/* -------------------------------------------------------------------- */
/*      Extract the data value in a way that depends on the number      */
/*      of bits in it.                                                  */
/* -------------------------------------------------------------------- */
        if( nNumBits == 0 )
        {
            nDataValue = 0;
        }
        else if( nNumBits == 1 )
        {
            nDataValue =
                (pabyValues[nValueBitOffset>>3] >> (nValueBitOffset&7)) & 0x1;
            nValueBitOffset++;
        }
        else if( nNumBits == 2 )
        {
            nDataValue =
                (pabyValues[nValueBitOffset>>3] >> (nValueBitOffset&7)) & 0x3;
            nValueBitOffset += 2;
        }
        else if( nNumBits == 4 )
        {
            nDataValue =
                (pabyValues[nValueBitOffset>>3] >> (nValueBitOffset&7)) & 0xf;
            nValueBitOffset += 4;
        }
        else if( nNumBits == 8 )
        {
            nDataValue = *pabyValues;
            pabyValues++;
        }
        else if( nNumBits == 16 )
        {
            nDataValue = 256 * *(pabyValues++);
            nDataValue += *(pabyValues++);
        }
        else if( nNumBits == 32 )
        {
            nDataValue = 256 * 256 * 256 * *(pabyValues++);
            nDataValue += 256 * 256 * *(pabyValues++);
            nDataValue += 256 * *(pabyValues++);
            nDataValue += *(pabyValues++);
        }
        else
        {
            printf( "nNumBits = %d\n", nNumBits );
            CPLAssert( FALSE );
            nDataValue = 0;
        }

/* -------------------------------------------------------------------- */
/*      Offset by the minimum value.                                    */
/* -------------------------------------------------------------------- */
        nDataValue += nDataMin;

/* -------------------------------------------------------------------- */
/*      Now apply to the output buffer in a type specific way.          */
/* -------------------------------------------------------------------- */
        if( nPixelsOutput + nRepeatCount > nMaxPixels )
        {
            CPLAssert( FALSE );
            nRepeatCount = nMaxPixels - nPixelsOutput;
        }
        
        if( nDataType == EPT_u8 )
        {
            int		i;
            
            for( i = 0; i < nRepeatCount; i++ )
            {
                CPLAssert( nDataValue < 256 );
                ((GByte *) pabyDest)[nPixelsOutput++] = (GByte)nDataValue;
            }
        }
        else if( nDataType == EPT_u16 )
        {
            int		i;
            
            for( i = 0; i < nRepeatCount; i++ )
            {
                ((GUInt16 *) pabyDest)[nPixelsOutput++] = (GUInt16)nDataValue;
            }
        }
        else if( nDataType == EPT_s16 )
        {
            int		i;
            
            for( i = 0; i < nRepeatCount; i++ )
            {
                ((GInt16 *) pabyDest)[nPixelsOutput++] = (GInt16)nDataValue;
            }
        }
        else if( nDataType == EPT_u32 )
        {
            int		i;
            
            for( i = 0; i < nRepeatCount; i++ )
            {
                ((GUInt32 *) pabyDest)[nPixelsOutput++] = (GUInt32)nDataValue;
            }
        }
        else if( nDataType == EPT_s32 )
        {
            int		i;
            
            for( i = 0; i < nRepeatCount; i++ )
            {
                ((GInt32 *) pabyDest)[nPixelsOutput++] = (GInt32)nDataValue;
            }
        }
        else if( nDataType == EPT_f32 )
        {
            int		i;
            float fDataValue;

            memcpy( &fDataValue, &nDataValue, 4);
            for( i = 0; i < nRepeatCount; i++ )
            {
                ((float *) pabyDest)[nPixelsOutput++] = fDataValue;
            }
        }
        else if( nDataType == EPT_u1 )
        {
            int		i;

            CPLAssert( nDataValue == 0 || nDataValue == 1 );
            
            if( nDataValue == 1 )
            {
                for( i = 0; i < nRepeatCount; i++ )
                {
                    pabyDest[nPixelsOutput>>3] |= (1 << (nPixelsOutput & 0x7));
                    nPixelsOutput++;
                }
            }
            else
            {
                for( i = 0; i < nRepeatCount; i++ )
                {
                    pabyDest[nPixelsOutput>>3] &= ~(1<<(nPixelsOutput & 0x7));
                    nPixelsOutput++;
                }
            }
        }
        else if( nDataType == EPT_u4 )
        {
            int		i;

            CPLAssert( nDataValue >= 0 && nDataValue < 16 );
            
            for( i = 0; i < nRepeatCount; i++ )
            {
                if( (nPixelsOutput & 0x1) == 0 )
                    pabyDest[nPixelsOutput>>1] = (GByte) nDataValue;
                else
                    pabyDest[nPixelsOutput>>1] |= (GByte) (nDataValue<<4);

                nPixelsOutput++;
            }
        }
        else
        {
            CPLError( CE_Failure, CPLE_AppDefined, 
                      "Attempt to uncompress an unsupported pixel data type.");
            return CE_Failure;
        }
    }

    return CE_None;
}