Browse Source

Puts quantization in another namespace and adds non-linear quantization out of the box.

master
Joshua Moerman 10 years ago
parent
commit
7d6ffce543
  1. 39
      wavelet/jcmp.cpp
  2. 15
      wavelet/jcmp_image.hpp
  3. 124
      wavelet/jcmp_quantization.hpp

39
wavelet/jcmp.cpp

@ -7,40 +7,29 @@
#include "jcmp.hpp"
#include "wavelet.hpp"
using namespace wvlt::V2;
// Can be any functions from R+ -> R
// as long as backward is its inverse
static double forward(double x){
//return std::log(1 + std::log(1 + x));
return std::log(x);
}
static double backward(double x){
//return std::exp(std::exp(x) - 1) - 1;
return std::exp(x);
}
using namespace wvlt;
using namespace jcmp;
// note: we take a copy, because we will modify it in place
static jcmp::image compress(std::vector<double> image, jcmp::uint width, double threshold, unsigned int& nzeros){
jcmp::uint height = image.size() / width;
static image compress(std::vector<double> img, uint width, double threshold, unsigned int& nzeros){
uint height = img.size() / width;
assert(is_pow_of_two(width));
assert(is_pow_of_two(height));
// wavelet transform in x-direction
for(unsigned int i = 0; i < height; ++i){
wavelet(&image[i*width], width, 1);
wavelet(&img[i*width], width, 1);
}
// wavelet transform in y-direction
for(unsigned int i = 0; i < width; ++i){
wavelet(&image[i], height, width);
wavelet(&img[i], height, width);
}
double min_abs = 10000.0;
double max_abs = 0.0;
for(auto& el : image){
for(auto& el : img){
auto absel = std::abs(el);
if(absel > threshold) {
min_abs = std::min(min_abs, absel);
@ -50,28 +39,28 @@ static jcmp::image compress(std::vector<double> image, jcmp::uint width, double
}
}
jcmp::quantization q(&forward, &backward, max_abs, min_abs);
auto q = quantization::get(quantization::type::logarithmic, max_abs, min_abs);
// save the principal coefficients
std::vector<jcmp::coefficient> v;
for(unsigned int y = 0; y < height; ++y){
for(unsigned int x = 0; x < width; ++x){
auto&& el = image[x + width*y];
if(el != 0) v.push_back({q.forwards(el), jcmp::uint(x), jcmp::uint(y)});
auto&& el = img[x + width*y];
if(el != 0) v.push_back({q.forwards(el), uint(x), uint(y)});
}
}
nzeros = v.size();
return jcmp::image(width, height, q.p, std::move(v));
return image(width, height, q.p, std::move(v));
}
static std::vector<double> decompress(jcmp::image in){
static std::vector<double> decompress(image in){
auto width = in.header.width;
auto height = in.header.height;
assert(is_pow_of_two(width));
assert(is_pow_of_two(height));
auto q = in.header.get_quantization(&forward, &backward);
auto q = in.header.get_quantization(quantization::type::logarithmic);
std::vector<double> image(width * height, 0.0);
@ -120,7 +109,7 @@ int main(){
// compress and decompress to see how we lost information
unsigned int nzeros = 0;
auto compressed_vec = decompress(compress(image_vec, width, 0.1, nzeros));
auto compressed_vec = decompress(compress(image_vec, width, 0.5, nzeros));
// output some information
std::cout << field("raw") << human_string(sizeof(uint8_t) * image_vec.size(), "b") << std::endl;

15
wavelet/jcmp_image.hpp

@ -13,10 +13,10 @@ namespace jcmp {
uint width;
uint height;
uint length;
quantize_params qp;
quantization::parameters qp;
header() = default;
header(uint width_, uint height_, uint length_, quantize_params const & p)
header(uint width_, uint height_, uint length_, quantization::parameters const & p)
: signature{'J', 'C', 'M', 'P'}
, width(width_)
, height(height_)
@ -24,13 +24,8 @@ namespace jcmp {
, qp(p)
{}
template <typename F>
quantization get_quantization(F const & f, F const & f_inv){
quantization q;
q.p = qp;
q.f = f;
q.f_inv = f_inv;
return q;
quantization::base get_quantization(quantization::type t){
return quantization::get(t, qp.f_max_abs, qp.f_min_abs, false);
}
};
@ -46,7 +41,7 @@ namespace jcmp {
image() = default;
image(uint width, uint height, quantize_params const & p, std::vector<coefficient> const & data_in)
image(uint width, uint height, quantization::parameters const & p, std::vector<coefficient> const & data_in)
: header(width, height, data_in.size(), p)
, data(data_in)
{}

124
wavelet/jcmp_quantization.hpp

@ -1,49 +1,91 @@
#pragma once
namespace jcmp {
// input in [0, 1], output [0, 255]
inline uint8_t clamp_round(double x){
if(x <= 0.0) return 0;
if(x >= 1.0) return 255;
return uint8_t(std::round(255.0*x));
}
#include <functional>
// parameters which need to be stored in the file
// to be able to dequantize
struct quantize_params {
double f_max_abs;
double f_min_abs;
};
// Quantization may be non linear (given by f and f_inv)
struct quantization {
std::function<double(double)> f;
std::function<double(double)> f_inv;
quantize_params p;
quantization() = default;
template <typename F>
quantization(F const & f_, F const & f_inv_, double max_abs, double min_abs)
: f(f_)
, f_inv(f_inv_)
, p{f(max_abs), f(min_abs)}
{}
uint8_t forwards(double x){
double y = std::abs(x);
y = (f(y) - p.f_min_abs) / (p.f_max_abs - p.f_min_abs);
if(x < 0) y = -y;
return clamp_round(0.5 * (y + 1.0));
namespace jcmp {
// Will quantize a double such that it fits in a uint8_t
// Also support non-linear quantization (which is usefull
// when the data is non-uniformly distributed)
namespace quantization {
// input in [0, 1], output [0, 255]
inline uint8_t clamp_round(double x){
if(x <= 0.0) return 0;
if(x >= 1.0) return 255;
return uint8_t(std::round(255.0*x));
}
double backwards(uint8_t x){
double y = 2.0 * x / 255.0 - 1.0;
y = (p.f_max_abs - p.f_min_abs) * y;
if(y < 0) return -f_inv(-y + p.f_min_abs);
return f_inv(y + p.f_min_abs);
// parameters which need to be stored in the file
// to be able to dequantize
// TODO: add type of quantization
struct parameters {
double f_max_abs;
double f_min_abs;
};
// the basic struct doing the quantization
struct base {
// f needs to be of type R+ -> R
// and F-inv should be its inverse
std::function<double(double)> f;
std::function<double(double)> f_inv;
parameters p;
uint8_t forwards(double x){
double y = std::abs(x);
y = (f(y) - p.f_min_abs) / (p.f_max_abs - p.f_min_abs);
if(x < 0) y = -y;
return clamp_round(0.5 * (y + 1.0));
}
double backwards(uint8_t x){
double y = 2.0 * x / 255.0 - 1.0;
y = (p.f_max_abs - p.f_min_abs) * y;
if(y < 0) return -f_inv(-y + p.f_min_abs);
return f_inv(y + p.f_min_abs);
}
};
// kinds of quantization which come out of the
// box with get(...)
enum class type {
linear,
logarithmic,
square_root
};
// non-overloaded versions
double id(double x) { return x; }
double log(double x) { return std::log(x); }
double exp(double x) { return std::exp(x); }
double sqrt(double x) { return std::sqrt(x); }
double sqr(double x) { return x*x; }
base get(type t, double max_abs, double min_abs, bool apply = true){
base b;
switch (t) {
case type::linear:
b.f = &id;
b.f_inv = &id;
break;
case type::logarithmic:
b.f = &log;
b.f_inv = &exp;
break;
case type::square_root:
b.f = &sqrt;
b.f_inv = &sqr;
break;
}
if(apply){
b.p.f_max_abs = b.f(max_abs);
b.p.f_min_abs = b.f(min_abs);
} else {
b.p.f_max_abs = max_abs;
b.p.f_min_abs = min_abs;
}
return b;
}
};
static_assert(sizeof(quantize_params) == 16, "struct not propery packed");
static_assert(sizeof(parameters) == 16, "struct not propery packed");
}
}