Image compression using wavelet
Wavelet compression is the generic name for a class of image coding methods that use a two-dimensional wavelet decomposition of an image to be encoded. Usually implied lossy compression. The article will not provide complex mathematical formulas, the whole theory can be read on the links at the bottom of the article. Here is only practice!
The JPEG algorithm, unlike the wavelet one, individually compresses each block of the original image of 8 by 8 pixels. As a result, a block structure may be noticeable at high degrees of compression on the reconstructed image. With wavelet compression, this problem does not arise, but other types of distortion may appear, having the appearance of “ghostly” ripples near sharp boundaries.
It is believed that such artifacts on average are less noticeable to the observer than the "squares" created by JPEG.
For example, we strongly compress the same image to approximately the same size:
')
At the beginning using JPEG:
(7959 bytes)
then the JPEG 2000 wavelet compression algorithm:
(7813 bytes)
In these examples, you can immediately see the nature of the introduced distortions by both algorithms. The image compressed using wavelet looks clearly better.
Another example of compressing another image , with a reduction in file size of up to 3kb each:
(Jpeg)
(JPEG 2000)
The method I implemented is in essence an analogue of the well-known JPEG 2000 wavelet compression algorithm, which has not replaced the popular JPEG. A few years ago I re-branded it into VisualBasic6.0, but especially for understanding the source code for most of the audience, I rewrote it in C #.
In the presented implementation, I used a fast lifting of a discrete biorthogonal CDF 9/7 wavelet used in the JPEG 2000 image compression algorithm. You can download the C implementation here .
In order not to inflate the code, I did not tolerate the ability to compress images of any size and with any aspect ratio. Therefore, it can compress only square images, the size of the sides 256x256, 512x512, 1024x1024, 2048x2048, etc. I also removed the memory optimization and saving / reading of the header of the file with the coefficients used in compression. Otherwise, the code would again be unnecessarily bloated.
I am not an ace in C #, so maybe there are places in the code where you could use the standard .Net methods.
Naturally the decompressor does everything in reverse order.
As I already wrote above, for the sake of readability of the code, I did not begin to save to the header file.
Therefore, in the run () method, the size of the decoded image and the coefficients for decoding the wavelet are rigidly clogged.
WaveleteCompression.zip
mikhanoid : 1.149604398, -0.4435068522, -0.8829110762, etc. - coefficients of the values of the functions of the wavelet basis at certain points. Read in general: Lifting scheme
JPEG difference
The JPEG algorithm, unlike the wavelet one, individually compresses each block of the original image of 8 by 8 pixels. As a result, a block structure may be noticeable at high degrees of compression on the reconstructed image. With wavelet compression, this problem does not arise, but other types of distortion may appear, having the appearance of “ghostly” ripples near sharp boundaries.
It is believed that such artifacts on average are less noticeable to the observer than the "squares" created by JPEG.
Example
For example, we strongly compress the same image to approximately the same size:
')
At the beginning using JPEG:
(7959 bytes)
then the JPEG 2000 wavelet compression algorithm:
(7813 bytes)
In these examples, you can immediately see the nature of the introduced distortions by both algorithms. The image compressed using wavelet looks clearly better.
Another example of compressing another image , with a reduction in file size of up to 3kb each:
(Jpeg)
(JPEG 2000)
What is it?
The method I implemented is in essence an analogue of the well-known JPEG 2000 wavelet compression algorithm, which has not replaced the popular JPEG. A few years ago I re-branded it into VisualBasic6.0, but especially for understanding the source code for most of the audience, I rewrote it in C #.
In the presented implementation, I used a fast lifting of a discrete biorthogonal CDF 9/7 wavelet used in the JPEG 2000 image compression algorithm. You can download the C implementation here .
The limitations of the code presented here are:
In order not to inflate the code, I did not tolerate the ability to compress images of any size and with any aspect ratio. Therefore, it can compress only square images, the size of the sides 256x256, 512x512, 1024x1024, 2048x2048, etc. I also removed the memory optimization and saving / reading of the header of the file with the coefficients used in compression. Otherwise, the code would again be unnecessarily bloated.
Main compression script:
- Load image from file;
- Convert the loaded image into a byte array of RGB values;
- Recode RGB to YCrCb with quantization of the resulting color components;
- We apply wavelet;
- We translate a multidimensional array into one-dimensional (flat);
- We compress the received stream by any available means (for example GZip).
Compressor Code (C #)
I am not an ace in C #, so maybe there are places in the code where you could use the standard .Net methods.
using System; using System.Collections.Generic; using System.Text; using System.Drawing; using System.Drawing.Imaging; using System.Runtime.InteropServices; using System.IO.Compression; using System.IO; namespace WaveleteCompression { class wvCompress { // public const int WV_LEFT_TO_RIGHT = 0; public const int WV_TOP_TO_BOTTOM = 1; public byte[] run(string path) { // Bitmap bmp = new Bitmap(path, true); // byte[, ,] b = this.BmpToBytes_Unsafe(bmp); // byte[] o = this.Compress(b, bmp.Width, bmp.Height); // RAW - FileStream f = new System.IO.FileStream(path + ".raw", FileMode.Create, FileAccess.Write); f.Write(o, 0, o.Length); f.Close(); // Gzip- // GZIP, 2 string outGZ = path + ".gz"; FileStream outfile = new FileStream(outGZ, FileMode.Create); GZipStream compressedzipStream = new GZipStream(outfile, CompressionMode.Compress, true); compressedzipStream.Write(o, 0, o.Length); compressedzipStream.Close(); // GZip- return o; } private byte[] Compress(byte[, ,] rgb, int cW, int cH) { // , int[] dwDiv = { 48, 32, 16, 16, 24, 24, 1, 1 }; int[] dwTop = { 24, 32, 24, 24, 24, 24, 32, 32 }; int SamplerDiv = 2, SamplerTop = 2; // Y, cr, cb int YPerec = 100, crPerec = 85, cbPerec = 85; int WVCount = 6; // // RGB YCrCb double[, ,] YCrCb = YCrCbEncode(rgb, cW, cH, YPerec, crPerec, cbPerec, cW, cH); // for (int z = 0; z < 3; z++) { // for (int dWave = 0; dWave < WVCount; dWave++) { int waveW = Convert.ToInt32(cW / Math.Pow(2, dWave)); int waveH = Convert.ToInt32(cH / Math.Pow(2, dWave)); if (z == 2) { // Y , // .. ( ), YCrCb = WaveletePack(YCrCb, z, waveW, waveH, dwDiv[dWave], dwTop[dWave], dWave); } else { YCrCb = WaveletePack(YCrCb, z, waveW, waveH, dwDiv[dWave] * SamplerDiv, dwTop[dWave] * SamplerTop, dWave); } } } // byte[] flattened = doPack(YCrCb, cW, cH, WVCount); return flattened; } /* Double Byte . Double Short. , , 255 */ private byte[] doPack(double[, ,] ImgData, int cW, int cH, int wDepth) { short Value; int lPos = 0; int size = cW * cH * 3; // short int intCount = 0; short[] shorts = new short[size]; byte[] Ret = new byte[size]; // - for(int d = wDepth-1; d >= 0; d--) { int wSize = (int)Math.Pow(2f, Convert.ToDouble(d)); int W = cW / wSize; int H = cH / wSize; int w2 = W / 2; int h2 = H / 2; // if (d == wDepth - 1) { for (int z = 0; z < 3; z++) { for (int j = 0; j < h2; j++) { for (int i = 0; i < w2; i++) { Value = (short)Math.Round(ImgData[z, i, j]); if ((Value >= -127) && (Value <= 127)) { Ret[lPos++] = Convert.ToByte(Value + 127); } else { Ret[lPos++] = 255; shorts[intCount++] = Value; } } } } } // + for (int z = 0; z < 3; z++) { for (int j = 0; j < H; j++) { for (int i = w2; i < W; i++) { Value = (short)Math.Round(ImgData[z, i, j]); if ((Value >= -127) && (Value <= 127)) { Ret[lPos++] = Convert.ToByte(Value + 127); } else { Ret[lPos++] = 255; shorts[intCount++] = Value; } } } } // for (int z = 0; z < 3; z++) { for (int j = h2; j < H; j++) { for (int i = 0; i < w2; i++) { Value = (short)Math.Round(ImgData[z, i, j]); if ((Value >= -127) && (Value <= 127)) { Ret[lPos++] = Convert.ToByte(Value + 127); } else { Ret[lPos++] = 255; shorts[intCount++] = Value; } } } } } // (byte[] short[]) int shortArraySize = intCount * 2; Array.Resize(ref Ret, Ret.Length + shortArraySize); Buffer.BlockCopy(shorts, 0, Ret, Ret.Length - shortArraySize, shortArraySize); // return Ret; } private double[, ,] WaveletePack(double[, ,] ImgArray, int Component, int cW, int cH, int dwDevider, int dwTop, int dwStep) { short Value; int cw2 = cW / 2; int cH2 = cH / 2; // double dbDiv = 1f / dwDevider; ImgArray = Wv(ImgArray, cW, cH, Component, WV_TOP_TO_BOTTOM); ImgArray = Wv(ImgArray, cH, cW, Component, WV_LEFT_TO_RIGHT); // for (int j = 0; j < cH; j++) { for (int i = 0; i < cW; i++) { if ((i >= cw2) || (j >= cH2)) { Value = (short)Math.Round(ImgArray[Component, i, j]); if (Value != 0) { int value2 = Value; if (value2 < 0) { value2 = -value2; } if (value2 < dwTop) { ImgArray[Component, i, j] = 0; } else { ImgArray[Component, i, j] = Value * dbDiv; } } } } } return ImgArray; } // CDF 9/7 private double[, ,] Wv(double[, ,] ImgArray, int n, int dwCh, int Component, int Side) { double a; int i, j, n2 = n / 2; double[] xWavelet = new double[n]; double[] tempbank = new double[n]; for (int dwPos = 0; dwPos < dwCh; dwPos++) { if (Side == WV_LEFT_TO_RIGHT) { for (j = 0; j < n; j++) { xWavelet[j] = ImgArray[Component, dwPos, j]; } } else if (Side == WV_TOP_TO_BOTTOM) { for (i = 0; i < n; i++) { xWavelet[i] = ImgArray[Component, i, dwPos]; } } // Predict 1 a = -1.586134342f; for (i = 1; i < n - 1; i += 2) { xWavelet[i] += a * (xWavelet[i - 1] + xWavelet[i + 1]); } xWavelet[n - 1] += 2 * a * xWavelet[n - 2]; // Update 1 a = -0.05298011854f; for (i = 2; i < n; i += 2) { xWavelet[i] += a * (xWavelet[i - 1] + xWavelet[i + 1]); } xWavelet[0] += 2 * a * xWavelet[1]; // Predict 2 a = 0.8829110762f; for (i = 1; i < n - 1; i += 2) { xWavelet[i] += a * (xWavelet[i - 1] + xWavelet[i + 1]); } xWavelet[n - 1] += 2 * a * xWavelet[n - 2]; // Update 2 a = 0.4435068522f; for (i = 2; i < n; i += 2) { xWavelet[i] += a * (xWavelet[i - 1] + xWavelet[i + 1]); } xWavelet[0] += 2 * a * xWavelet[1]; // Scale a = 1f / 1.149604398f; j = 0; // "" // "" if (Side == WV_LEFT_TO_RIGHT) { for (i = 0; i < n2; i++) { ImgArray[Component, dwPos, i] = xWavelet[j++] / a; ImgArray[Component, dwPos, n2 + i] = xWavelet[j++] * a; } } else if (Side == WV_TOP_TO_BOTTOM) { for (i = 0; i < n2; i++) { ImgArray[Component, i, dwPos] = xWavelet[j++] / a; ImgArray[Component, n2 + i, dwPos] = xWavelet[j++] * a; } } } return ImgArray; } // RGB YCrCb private double[, ,] YCrCbEncode(byte[, ,] BytesRGB, int cW, int cH, double Ydiv, double Udiv, double Vdiv, int oW, int oH) { double vr, vg, vb; double kr = 0.299, kg = 0.587, kb = 0.114, kr1 = -0.1687, kg1 = 0.3313, kb1 = 0.5, kr2 = 0.5, kg2 = 0.4187, kb2 = 0.0813; Ydiv = Ydiv / 100f; Udiv = Udiv / 100f; Vdiv = Vdiv / 100f; double[, ,] YCrCb = new double[3, cW, cH]; for (int j = 0; j < oH; j++) { for (int i = 0; i < oW; i++) { vb = (double)BytesRGB[0, i, j]; vg = (double)BytesRGB[1, i, j]; vr = (double)BytesRGB[2, i, j]; YCrCb[2, i, j] = (kr * vr + kg * vg + kb * vb) * Ydiv; YCrCb[1, i, j] = (kr1 * vr - kg1 * vg + kb1 * vb + 128) * Udiv; YCrCb[0, i, j] = (kr2 * vr - kg2 * vg - kb2 * vb + 128) * Udiv; } } return YCrCb; } private unsafe byte[, ,] BmpToBytes_Unsafe(Bitmap bmp) { BitmapData bData = bmp.LockBits(new Rectangle(new Point(), bmp.Size), ImageLockMode.ReadOnly, PixelFormat.Format24bppRgb); // number of bytes in the bitmap int byteCount = bData.Stride * bmp.Height; byte[] bmpBytes = new byte[byteCount]; Marshal.Copy(bData.Scan0, bmpBytes, 0, byteCount); // Copy the locked bytes from memory // don't forget to unlock the bitmap!! bmp.UnlockBits(bData); byte[, ,] ret = new byte[3, bmp.Width, bmp.Height]; for (int z = 0; z < 3; z++) { for (int i = 0; i < bmp.Width; i++) { for (int j = 0; j < bmp.Height; j++) { ret[z, i, j] = bmpBytes[j * bmp.Width * 3 + i * 3 + z]; } } } return ret; } } } Decompressor code
Naturally the decompressor does everything in reverse order.
As I already wrote above, for the sake of readability of the code, I did not begin to save to the header file.
Therefore, in the run () method, the size of the decoded image and the coefficients for decoding the wavelet are rigidly clogged.
using System; using System.Collections.Generic; using System.Text; namespace WaveleteCompression { class wvDecompress { // public const int WV_LEFT_TO_RIGHT = 0; public const int WV_TOP_TO_BOTTOM = 1; public byte[] run(byte[] compressed) { int z; int dwDepth = 6; // ( , ) // int w = 512; int h = 512; // (header) int[] dwDiv = { 48, 32, 16, 16, 24, 24, 1, 1 }, dwTop = { 24, 32, 24, 24, 24, 24, 32, 32 }; int SamplerDiv = 2, YPerec = 100, crPerec = 85, cbPerec = 85; double[,,] yuv = doUnPack(compressed, w, h, dwDepth); // for(z = 0; z < 2; z++) { for(int dWave = dwDepth - 1; dWave >= 0; dWave--) { int w2 = Convert.ToInt32(w / Math.Pow(2, dWave)); int h2 = Convert.ToInt32(h / Math.Pow(2, dWave)); WaveleteUnPack(yuv, z, w2, h2, dwDiv[dWave] * SamplerDiv); } } z = 2; for(int dWave = dwDepth - 1; dWave >= 0; dWave--) { int w2 = Convert.ToInt32(w / Math.Pow(2, dWave)); int h2 = Convert.ToInt32(h / Math.Pow(2, dWave)); WaveleteUnPack(yuv, z, w2, h2, dwDiv[dWave]); } // YCrCb + byte[] rgb_flatened = this.YCrCbDecode(yuv, w, h, YPerec, crPerec, cbPerec); return rgb_flatened; } // DoPack wvCompress. // (short)double- byte[] private static double[, ,] doUnPack(byte[] Bytes, int cW, int cH, int dwDepth) { int lPos = 0; byte Value; int intIndex = 0; // int size = cW * cH * 3; // double[,,] ImgData = new double[3, cW, cH]; int shortsLength = Bytes.Length - size; short[] shorts = new short[shortsLength / 2]; Buffer.BlockCopy(Bytes, size, shorts, 0, shortsLength); for (int d = dwDepth - 1; d >= 0; d--) { int wSize = (int)Math.Pow(2, d); int W = cW / wSize; int H = cH / wSize; int w2 = W / 2; int h2 = H / 2; // if (d == dwDepth - 1) { for (int z = 0; z < 3; z++) { for (int j = 0; j < h2; j++) { for (int i = 0; i < w2; i++) { Value = Bytes[lPos++]; if(Value == 255) { ImgData[z, i, j] = shorts[intIndex++]; } else { ImgData[z, i, j] = Value - 127; } } } } } // + for (int z = 0; z < 3; z++) { for (int j = 0; j < H; j++) { for (int i = w2; i < W; i++) { Value = Bytes[lPos++]; if(Value == 255) { ImgData[z, i, j] = shorts[intIndex++]; } else { ImgData[z, i, j] = Value - 127; } } } } // for (int z = 0; z < 3; z++) { for (int j = h2; j < H; j++) { for (int i = 0; i < w2; i++) { Value = Bytes[lPos++]; if (Value == 255) { ImgData[z, i, j] = shorts[intIndex++]; } else { ImgData[z, i, j] = Value - 127; } } } } } // return ImgData; } // private void WaveleteUnPack(double[,,] ImgArray, int Component, int cW, int cH, int dwDevider) { int cw2 = cW / 2, ch2 = cH / 2; double dbDiv = 1f / dwDevider; // for(int i = 0; i < cW; i++) { for(int j = 0; j < cH; j++) { if ((i >= cw2) || (j >= ch2)) { if (ImgArray[Component, i, j] != 0) { ImgArray[Component, i, j] /= dbDiv; } } } } // for(int i = 0; i < cW; i++) { reWv(ref ImgArray, cH, Component, i, WV_LEFT_TO_RIGHT); } for(int j = 0; j < cH; j++) { reWv(ref ImgArray, cW, Component, j, WV_TOP_TO_BOTTOM); } } // CDF 9/7 private void reWv(ref double[, ,] shorts, int n, int z, int dwPos, int Side) { double a; double[] xWavelet = new double[n]; double[] tempbank = new double[n]; if(Side == WV_LEFT_TO_RIGHT) { for(int j = 0; j < n; j++) { xWavelet[j] = shorts[z, dwPos, j]; } } else if (Side == WV_TOP_TO_BOTTOM) { for(int i = 0; i < n; i++) { xWavelet[i] = shorts[z, i, dwPos]; } } for(int i = 0; i < n / 2; i++) { tempbank[i * 2] = xWavelet[i]; tempbank[i * 2 + 1] = xWavelet[i + n / 2]; } for(int i = 0; i < n; i++) { xWavelet[i] = tempbank[i]; } // Undo scale a = 1.149604398f; for(int i = 0; i < n; i++) { if(i % 2 != 0) { xWavelet[i] = xWavelet[i] * a; } else { xWavelet[i] = xWavelet[i] / a; } } // Undo update 2 a = -0.4435068522f; for (int i = 2; i < n; i += 2) { xWavelet[i] = xWavelet[i] + a * (xWavelet[i - 1] + xWavelet[i + 1]); } xWavelet[0] = xWavelet[0] + 2 * a * xWavelet[1]; // Undo predict 2 a = -0.8829110762f; for (int i = 1; i < n - 1; i += 2) { xWavelet[i] = xWavelet[i] + a * (xWavelet[i - 1] + xWavelet[i + 1]); } xWavelet[n - 1] = xWavelet[n - 1] + 2 * a * xWavelet[n - 2]; // Undo update 1 a = 0.05298011854f; for (int i = 2; i < n; i += 2) { xWavelet[i] = xWavelet[i] + a * (xWavelet[i - 1] + xWavelet[i + 1]); } xWavelet[0] = xWavelet[0] + 2 * a * xWavelet[1]; // Undo predict 1 a = 1.586134342f; for (int i = 1; i < n - 1; i += 2) { xWavelet[i] = xWavelet[i] + a * (xWavelet[i - 1] + xWavelet[i + 1]); } xWavelet[n - 1] = xWavelet[n - 1] + 2 * a * xWavelet[n - 2]; if(Side == WV_LEFT_TO_RIGHT) { for (int j = 0; j < n; j++) { shorts[z, dwPos, j] = xWavelet[j]; } } else if(Side == WV_TOP_TO_BOTTOM) { for(int i = 0; i < n; i++) { shorts[z, i, dwPos] = xWavelet[i]; } } } // YCrCb RGB private byte[] YCrCbDecode(double[, ,] yuv, int w, int h, double Ydiv, double Udiv, double Vdiv) { byte[] bytes_flat = new byte[3 * w * h]; double vr, vg, vb; double vY, vCb, vCr; Ydiv = Ydiv / 100f; Udiv = Udiv / 100f; Vdiv = Vdiv / 100f; for(int j = 0; j < h; j++) { for (int i = 0; i < w ; i++) { vCr = yuv[0, i, j] / Vdiv; vCb = yuv[1, i, j] / Udiv; vY = yuv[2, i, j] / Ydiv; vr = vY + 1.402f * (vCr - 128f); vg = vY - 0.34414f * (vCb - 128f) - 0.71414f * (vCr - 128f); vb = vY + 1.722f * (vCb - 128f); if (vr > 255) {vr = 255;} if (vg > 255) {vg = 255;} if (vb > 255) {vb = 255;} if (vr < 0) {vr = 0;} if (vg < 0) {vg = 0;} if (vb < 0) {vb = 0;} bytes_flat[j * w * 3 + i * 3 + 0] = (byte)vb; bytes_flat[j * w * 3 + i * 3 + 1] = (byte)vg; bytes_flat[j * w * 3 + i * 3 + 2] = (byte)vr; } } return bytes_flat; } } } C # project source codes
WaveleteCompression.zip
Links
UPD:mikhanoid : 1.149604398, -0.4435068522, -0.8829110762, etc. - coefficients of the values of the functions of the wavelet basis at certain points. Read in general: Lifting scheme
Source: https://habr.com/ru/post/128913/
All Articles