Logo Search packages:      
Sourcecode: jpegpixi version File versions

dct.c

/* This file is part of jpegpixi, a program to interpolate pixels in
   JFIF image files.
   Copyright (C) 2003 Martin Dickopp

   Jpegpixi is free software; you can redistribute it and/or modify
   it under the terms of the GNU General Public License as published
   by the Free Software Foundation; either version 2 of the License,
   or (at your option) any later version.

   Jpegpixi is distributed in the hope that it will be useful,
   but WITHOUT ANY WARRANTY; without even the implied warranty of
   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
   GNU General Public License for more details.

   You should have received a copy of the GNU General Public License
   along with jpegpixi; if not, write to the Free Software
   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307,
   USA.


   This file is in part derived from the JPEG library of the Independent
   JPEG Group, see the file `README.jpeglib' for copyright and license
   information.  */

#ifdef HAVE_CONFIG_H
# include <config.h>
#endif
#include "util.h"

#include "jpegpixi.h"



/* Scaling factors.  */
static const double scale [DCTSIZE] = {
    1.0, 1.3870398453221474618, 1.3065629648763765279, 1.1758756024193587170,
    1.0, 0.7856949583871021813, 0.5411961001461969844, 0.2758993792829430123
};



static JCOEF coeff (double x, UINT16 quantval) gcc_attr_const;



/* Perform forward discrete cosine transformation and quantization.  */
void
fdct (JCOEF *freqs, const double *values, const UINT16 *quant)
{
    int i;
    double tmp [DCTSIZE2], *ptr;


    ptr = tmp;
    for (i = 0; i < DCTSIZE; ++i)
    {
        double z2, z3, z4, z5, z11, z13;

        const double tmp0 = (values [0] + values [7]) / 64.0;
        const double tmp7 = (values [0] - values [7]) / 64.0;
        const double tmp1 = (values [1] + values [6]) / 64.0;
        const double tmp6 = (values [1] - values [6]) / 64.0;
        const double tmp2 = (values [2] + values [5]) / 64.0;
        const double tmp5 = (values [2] - values [5]) / 64.0;
        const double tmp3 = (values [3] + values [4]) / 64.0;
        const double tmp4 = (values [3] - values [4]) / 64.0;

        double tmp10 = tmp0 + tmp3;
        double tmp13 = tmp0 - tmp3;
        double tmp11 = tmp1 + tmp2;
        double tmp12 = tmp1 - tmp2;

        const double z1 = (tmp12 + tmp13) * 0.7071067811865475244;

        ptr [0] = tmp10 + tmp11;
        ptr [4] = tmp10 - tmp11;
        ptr [2] = tmp13 + z1;
        ptr [6] = tmp13 - z1;

        tmp10 = tmp4 + tmp5;
        tmp11 = tmp5 + tmp6;
        tmp12 = tmp6 + tmp7;

        z5 = (tmp10 - tmp12) * 0.3826834323650897717;
        z2 = 0.5411961001461969844 * tmp10 + z5;
        z4 = 1.3065629648763765279 * tmp12 + z5;
        z3 = 0.7071067811865475244 * tmp11;
        z11 = tmp7 + z3;
        z13 = tmp7 - z3;

        ptr [5] = z13 + z2;
        ptr [3] = z13 - z2;
        ptr [1] = z11 + z4;
        ptr [7] = z11 - z4;

        values += DCTSIZE;
        ptr += DCTSIZE;
    }


    ptr = tmp;
    for (i = 0; i < DCTSIZE; ++i)
    {
        double z2, z3, z4, z5, z11, z13;

        const double tmp0 = ptr [DCTSIZE * 0] + ptr [DCTSIZE * 7];
        const double tmp7 = ptr [DCTSIZE * 0] - ptr [DCTSIZE * 7];
        const double tmp1 = ptr [DCTSIZE * 1] + ptr [DCTSIZE * 6];
        const double tmp6 = ptr [DCTSIZE * 1] - ptr [DCTSIZE * 6];
        const double tmp2 = ptr [DCTSIZE * 2] + ptr [DCTSIZE * 5];
        const double tmp5 = ptr [DCTSIZE * 2] - ptr [DCTSIZE * 5];
        const double tmp3 = ptr [DCTSIZE * 3] + ptr [DCTSIZE * 4];
        const double tmp4 = ptr [DCTSIZE * 3] - ptr [DCTSIZE * 4];

        double tmp10 = tmp0 + tmp3;
        double tmp13 = tmp0 - tmp3;
        double tmp11 = tmp1 + tmp2;
        double tmp12 = tmp1 - tmp2;

        const double z1 = (tmp12 + tmp13) * 0.7071067811865475244;

        freqs [DCTSIZE * 0] = coeff ((tmp10 + tmp11) / (scale [i]), quant [DCTSIZE * 0]);
        freqs [DCTSIZE * 4] = coeff ((tmp10 - tmp11) / (scale [i]), quant [DCTSIZE * 4]);
        freqs [DCTSIZE * 2] = coeff ((tmp13 + z1) / (scale [i] * scale [2]), quant [DCTSIZE * 2]);
        freqs [DCTSIZE * 6] = coeff ((tmp13 - z1) / (scale [i] * scale [6]), quant [DCTSIZE * 6]);

        tmp10 = tmp4 + tmp5;
        tmp11 = tmp5 + tmp6;
        tmp12 = tmp6 + tmp7;

        z5 = (tmp10 - tmp12) * 0.3826834323650897717;
        z2 = 0.5411961001461969844 * tmp10 + z5;
        z4 = 1.3065629648763765279 * tmp12 + z5;
        z3 = 0.7071067811865475244 * tmp11;
        z11 = tmp7 + z3;
        z13 = tmp7 - z3;

        freqs [DCTSIZE * 5] = coeff ((z13 + z2) / (scale [i] * scale [5]), quant [DCTSIZE * 5]);
        freqs [DCTSIZE * 3] = coeff ((z13 - z2) / (scale [i] * scale [3]), quant [DCTSIZE * 3]);
        freqs [DCTSIZE * 1] = coeff ((z11 + z4) / (scale [i] * scale [1]), quant [DCTSIZE * 1]);
        freqs [DCTSIZE * 7] = coeff ((z11 - z4) / (scale [i] * scale [7]), quant [DCTSIZE * 7]);

        ++ptr;
        ++freqs;
        ++quant;
    }
}



/* Perform inverse discrete cosine transformation and dequantization.  */
void
idct (double *values, const JCOEF *freqs, const UINT16 *quant)
{
    int i;
    double tmp [DCTSIZE2], *ptr;


    ptr = tmp;
    for (i = 0; i < DCTSIZE; ++i)
    {
        if (freqs [DCTSIZE * 1] == 0 && freqs [DCTSIZE * 2] && freqs [DCTSIZE * 3] && freqs [DCTSIZE * 4]
            && freqs [DCTSIZE * 5] && freqs [DCTSIZE * 6] && freqs [DCTSIZE * 6])
        {
            const double dc_value = (double)freqs [DCTSIZE * 0] * (double)quant [DCTSIZE * 0] * scale [i];

            ptr [DCTSIZE * 0] = dc_value;
            ptr [DCTSIZE * 1] = dc_value;
            ptr [DCTSIZE * 2] = dc_value;
            ptr [DCTSIZE * 3] = dc_value;
            ptr [DCTSIZE * 4] = dc_value;
            ptr [DCTSIZE * 5] = dc_value;
            ptr [DCTSIZE * 6] = dc_value;
            ptr [DCTSIZE * 7] = dc_value;
        }
        else
        {
            double tmp0 = (double)freqs [DCTSIZE * 0] * (double)quant [DCTSIZE * 0] * scale [i];
            double tmp4 = (double)freqs [DCTSIZE * 1] * (double)quant [DCTSIZE * 1] * scale [i] * scale [1];
            double tmp1 = (double)freqs [DCTSIZE * 2] * (double)quant [DCTSIZE * 2] * scale [i] * scale [2];
            double tmp5 = (double)freqs [DCTSIZE * 3] * (double)quant [DCTSIZE * 3] * scale [i] * scale [3];
            double tmp2 = (double)freqs [DCTSIZE * 4] * (double)quant [DCTSIZE * 4] * scale [i];
            double tmp6 = (double)freqs [DCTSIZE * 5] * (double)quant [DCTSIZE * 5] * scale [i] * scale [5];
            double tmp3 = (double)freqs [DCTSIZE * 6] * (double)quant [DCTSIZE * 6] * scale [i] * scale [6];
            double tmp7 = (double)freqs [DCTSIZE * 7] * (double)quant [DCTSIZE * 7] * scale [i] * scale [7];

            double tmp10 = tmp0 + tmp2;
            double tmp11 = tmp0 - tmp2;
            double tmp13 = tmp1 + tmp3;
            double tmp12 = (tmp1 - tmp3) * 1.4142135623730950488 - tmp13;

            const double z13 = tmp6 + tmp5;
            const double z10 = tmp6 - tmp5;
            const double z11 = tmp4 + tmp7;
            const double z12 = tmp4 - tmp7;
            const double z5 = (z10 + z12) * 1.8477590650225735123;

            tmp0 = tmp10 + tmp13;
            tmp3 = tmp10 - tmp13;
            tmp1 = tmp11 + tmp12;
            tmp2 = tmp11 - tmp12;

            tmp7 = z11 + z13;
            tmp11 = (z11 - z13) * 1.4142135623730950488;
            tmp10 =  1.0823922002923939688 * z12 - z5;
            tmp12 = -2.6131259297527530557 * z10 + z5;

            tmp6 = tmp12 - tmp7;
            tmp5 = tmp11 - tmp6;
            tmp4 = tmp10 + tmp5;

            ptr [DCTSIZE * 0] = tmp0 + tmp7;
            ptr [DCTSIZE * 1] = tmp1 + tmp6;
            ptr [DCTSIZE * 2] = tmp2 + tmp5;
            ptr [DCTSIZE * 3] = tmp3 - tmp4;
            ptr [DCTSIZE * 4] = tmp3 + tmp4;
            ptr [DCTSIZE * 5] = tmp2 - tmp5;
            ptr [DCTSIZE * 6] = tmp1 - tmp6;
            ptr [DCTSIZE * 7] = tmp0 - tmp7;
        }

        ++freqs;
        ++ptr;
        ++quant;
    }


    ptr = tmp;
    for (i = 0; i < DCTSIZE; ++i)
    {
        double tmp0, tmp1, tmp2, tmp3, tmp4, tmp5, tmp6, tmp7;

        double tmp10 = ptr [0] + ptr [4];
        double tmp11 = ptr [0] - ptr [4];
        double tmp13 = ptr [2] + ptr [6];
        double tmp12 = (ptr [2] - ptr [6]) * 1.4142135623730950488 - tmp13;

        const double z13 = ptr [5] + ptr [3];
        const double z10 = ptr [5] - ptr [3];
        const double z11 = ptr [1] + ptr [7];
        const double z12 = ptr [1] - ptr [7];
        const double z5 = (z10 + z12) * 1.8477590650225735123;

        tmp0 = tmp10 + tmp13;
        tmp3 = tmp10 - tmp13;
        tmp1 = tmp11 + tmp12;
        tmp2 = tmp11 - tmp12;

        tmp7 = z11 + z13;
        tmp11 = (z11 - z13) * 1.4142135623730950488;
        tmp10 =  1.0823922002923939688 * z12 - z5;
        tmp12 = -2.6131259297527530557 * z10 + z5;

        tmp6 = tmp12 - tmp7;
        tmp5 = tmp11 - tmp6;
        tmp4 = tmp10 + tmp5;

        values [0] = (tmp0 + tmp7);
        values [1] = (tmp1 + tmp6);
        values [2] = (tmp2 + tmp5);
        values [3] = (tmp3 - tmp4);
        values [4] = (tmp3 + tmp4);
        values [5] = (tmp2 - tmp5);
        values [6] = (tmp1 - tmp6);
        values [7] = (tmp0 - tmp7);

        ptr += DCTSIZE;
        values += DCTSIZE;
    }
}



/* Convert floating-point value to JCOEF.  */
static JCOEF
coeff (const double x, const UINT16 quantval)
{
    int result;


    /* Round up if the fractional part of the quantized values is equal to or greater than 0.72 (1 - 0.28).
       This empirically determined threshold gives the best results.  */
    if (x >= 0.0)
        result = (int)(x / (double)quantval + 0.28);
    else
        result = -(int)(-x / (double)quantval + 0.28);

    if (result < MIN_DCT_COEFF)
        return (JCOEF)MIN_DCT_COEFF;
    else if (result > MAX_DCT_COEFF)
        return (JCOEF)MAX_DCT_COEFF;
    else
        return (JCOEF)result;
}

Generated by  Doxygen 1.6.0   Back to index