1
\$\begingroup\$

The aim of my code is to decode an image format that is based on the JPEG chain of compression/decompression, however it is not compatible with the default JPEG flow as far as I know, since all libraries I have tried fail to properly decode the data. I am only interested decompression in this case. It follows the standard pattern:

  • Read Huffman values -> Like normal JPEG
  • unzigzag -> Like normal JPEG
  • Dequantize -> Like normal JPEG
  • IDCT -> Almost like normal JPEG, but different range/clamping
  • Color Space conversion -> Custom, not YCbCr

For one 8x8 except for the last step that looks like this right now:

int16_t processBlock(int16_t prevDc, BitStream &stream, const tHuffTable &dcTable, const tHuffTable &acTable,
                     float *quantTable, bool isLuminance, int16_t *outBlock) {
    int16_t workBlock[64] = {0};
    int16_t curDc = decodeBlock(stream, workBlock, dcTable, acTable, prevDc); 
    unzigzag(workBlock);
    dequantize(workBlock, quantTable);
    idct(outBlock, workBlock, isLuminance);
    return curDc;
}

after this the outBlock is treated by the color space conversion based on the image type.

What I want to optimize is the overall performance. The entire image is decompressed in the following way with 4 luminance blocks for component 1, 1 chrominance block for component 2 and 1 chrominance block for component 3. There are 4 more blocks for another luminance component, but I dont know what it is used for, so we can ignore it. The code looks like this:

void decodeImageType0(
        uint32_t width,
        uint32_t height,
        std::vector<uint8_t> &outData,
        BitStream &stream,
        const tHuffTable &dcLumTable,
        const tHuffTable &acLumTable,
        const tHuffTable &dcCromTable,
        const tHuffTable &acCromTable,
        float *lumQuant[4],
        float *cromQuant[4]) {
    int16_t lum0[4][64]{};
    int16_t lum1[4][64]{};
    int16_t crom0[64]{};
    int16_t crom1[64]{};
    uint32_t colorBlock[16 * 16]{};

    const auto actualHeight = ((height + 15) / 16) * 16;
    const auto actualWidth = ((width + 15) / 16) * 16;

    int16_t prevDc[4] = {0};
    for (auto y = 0; y < (actualHeight / 16); ++y) {
        for (auto x = 0; x < (actualWidth / 16); ++x) {
            for (auto &lum : lum0) {
                prevDc[0] = processBlock(prevDc[0], stream, dcLumTable, acLumTable, lumQuant[0], true, lum);
            }
            prevDc[1] = processBlock(prevDc[1], stream, dcCromTable, acCromTable, cromQuant[1], false, crom0);
            prevDc[2] = processBlock(prevDc[2], stream, dcCromTable, acCromTable, cromQuant[2], false, crom1);
            for (auto &lum : lum1) {
                prevDc[3] = processBlock(prevDc[3], stream, dcLumTable, acLumTable, lumQuant[3], true, lum);
            }

            decodeColorBlockType0(lum0, lum1, crom0, crom1, colorBlock);
            for (auto row = 0; row < 16; ++row) {
                if(y * 16 + row >= height || x * 16 >= width) {
                    continue;
                }

                const auto numPixels = std::min(16u, width - x * 16);
                memcpy(outData.data() + (y * 16 + row) * width * 4 + x * 16 * 4, &colorBlock[row * 16], numPixels * 4);
            }
        }
    }
}

Now my measurements have shown that over 80% of the time is spent inside the idct function, so this is where I want to optimize. The function looks like this, after I applied what I could think of to optimize it. I have created a cache of the static coefficients used in the IDCT process which significantly improved performance, but I hope there is still room for more, for example nanojpg is 3 times faster (however with invalid results).

float idctHelper(const int16_t *inBlock, int32_t u, int32_t v, int32_t blockWidth, int32_t blockHeight) {
    glm::vec<4, float, glm::packed_lowp> vec3{};

    float result = 0.0f;
    for (auto y = 0; y < blockHeight; ++y) {
        for (auto x = 0; x < blockWidth; x += 4) {
            const auto idx = (v * 8 + u) * 64 + y * 8 + x;
            vec3 = glm::vec<4, float, glm::packed_lowp>(inBlock[y * blockWidth + x], inBlock[y * blockWidth + x + 1], inBlock[y * blockWidth + x + 2], inBlock[y * blockWidth + x + 3]) *
                    glm::vec<4, float, glm::packed_lowp>(idctLookup[idx], idctLookup[idx + 1], idctLookup[idx + 2], idctLookup[idx + 3]);
            result += vec3.x + vec3.y + vec3.z + vec3.w;
        }
    }

    return result;
}

template<typename T, typename U = T>
U clamp(T value, T min, T max) {
    return static_cast<U>(std::min<T>(std::max<T>(value, min), max));
}

void idct(int16_t *outBlock, int16_t *inBlock, bool isLuminance, int32_t blockWidth = 8, int32_t blockHeight = 8) {
    for (auto y = 0; y < blockHeight; ++y) {
        for (auto x = 0; x < blockWidth; ++x) {
            auto value = static_cast<int16_t>(std::round(
                    0.25f * idctHelper(inBlock, x, y, blockWidth, blockHeight)));
            if (isLuminance) {
                value = clamp<int16_t>(static_cast<int16_t>(value + 128), 0, 255);
            } else {
                value = clamp<int16_t>(value, -256, 255);
            }

            outBlock[y * blockWidth + x] = value;
        }
    }
}

This is the cache that is created once at the start of the application outside of time measurements:

float alphaFunction(int32_t n) {
    static float INV_SQRT_2 = 1.0f / sqrtf(2.0f);

    if (n == 0) {
        return INV_SQRT_2;
    } else {
        return 1;
    }
}
        for (auto u = 0; u < 8; ++u) {
            for (auto v = 0; v < 8; ++v) {
                for (auto x = 0; x < 8; ++x) {
                    for (auto y = 0; y < 8; ++y) {
                        idctLookup[(v * 8 + u) * 64 + y * 8 + x] = alphaFunction(x) * alphaFunction(y) *
                                                                   cosf((2 * u + 1) * x * (float) M_PI / 16.0f) *
                                                                   cosf((2 * v + 1) * y * (float) M_PI / 16.0f);
                    }
                }
            }
        }
\$\endgroup\$
  • 1
    \$\begingroup\$ The code at the end with the four nested for loops — is that part of a function or what? \$\endgroup\$ – 200_success Jan 13 at 4:17
  • \$\begingroup\$ I have edited the question. It is created only once at the start of the application and is not part of the measured times. \$\endgroup\$ – Cromon Jan 13 at 8:10
0
\$\begingroup\$

DCT can be implemented with a fast algorithm. I understand you implemented it with a matrix multiplication, which is much less efficient.

DCT algorithms are similar to FFT ones. You can find many references with a simple search, on Wikipedia for example

\$\endgroup\$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.