00001
00008 #ifndef MAXLLOYDQUANTISIERER_H_
00009 #define MAXLLOYDQUANTISIERER_H_
00010
00011 #include <cmath>
00012 #include <itpp/itbase.h>
00013 #include <itpp/stat/histogram.h>
00014
00015 using namespace std;
00016
00017 namespace itpp {
00018
00023 template<class T>class MaxLloydQuantizer {
00024 public:
00025
00029 MaxLloydQuantizer() { initialized = trained = false; }
00030
00035 MaxLloydQuantizer(unsigned int R) { setTargetRate(R); }
00036
00045 MaxLloydQuantizer(unsigned int R,itpp::Vec<T> trainingdata,double epsilon=1e-3,unsigned int max_iter=1e3,unsigned hist_size=1000) { setTargetRate(R); train(trainingdata,epsilon,max_iter,hist_size); }
00046
00050 virtual ~MaxLloydQuantizer() {}
00051
00059 bool train(const itpp::Vec<T> &trainingdata,double epsilon=1e-3,unsigned int max_iter=1e3);
00060
00061 bool train2(const itpp::Vec<T> &trainingdata,double epsilon=1e-3,unsigned int max_iter=1e3);
00062
00067 unsigned int getNumberOfQuantizerLevels() const { return N; }
00068
00073 unsigned int getTargetRate() const { return Rate; }
00074
00079 bool isTrained() const { return trained; }
00080
00085 bool isInitialized() const { return initialized; }
00086
00092 void setTargetRate(unsigned int R)
00093 {
00094 Rate=R;
00095 N=1<<R;
00096 d.set_length(N+1);
00097 r.set_length(N);
00098 trained=false;
00099 initialized=true;
00100 }
00101
00109 unsigned int quantize(T x);
00110
00116 T inv_quant(unsigned int level);
00117
00118 protected:
00119 itpp::Vec<T> d;
00120 itpp::Vec<T> r;
00121 unsigned int Rate;
00122 unsigned int N;
00123 bool trained;
00124 bool initialized;
00125 };
00126
00127
00128 template<class T> T MaxLloydQuantizer<T>::inv_quant(unsigned int level)
00129 {
00130 #ifdef DEBUG
00131 if(!initialized) cerr << "MaxLloydQuantizer<T>::inv_quant() called before initializing the quantizer!" << endl;
00132 if(!trained) cerr << "MaxLloydQuantizer<T>::inv_quant() called before training the quantizer!" << endl;
00133 if(level >= N) if(!trained) cerr << "MaxLloydQuantizer<T>::inv_quant(level): level out of range!" << endl;
00134 #endif
00135 return r(level);
00136 }
00137
00138
00139
00140
00141
00142
00143
00144
00145
00146
00147
00148
00149
00150
00151
00152
00153
00154
00155
00156
00157
00158
00159
00160
00161
00162
00163
00164
00165
00166
00167
00168
00169
00170
00171
00172
00173
00174
00175
00176
00177
00178
00179
00180
00181
00182
00183
00184
00185
00186
00187
00188
00189
00190
00191
00192
00193
00194
00195
00196
00197
00198
00199 template<class T> bool MaxLloydQuantizer<T>::train(const itpp::Vec<T> &trainingdata,double epsilon,unsigned int max_iter)
00200 {
00201 #ifdef DEBUG
00202 if(!initialized) cerr << "MaxLloydQuantizer<T>::train() called before initializing the quantizer!" << endl;
00203 #endif
00204 if(trainingdata.length()==0) return false;
00205 T minv,maxv;
00206 minv=maxv=trainingdata(0);
00207 for(const T* ptr=trainingdata._data();ptr < trainingdata._data()+trainingdata.length();ptr++)
00208 {
00209 register T val = *ptr;
00210 minv = ((val < minv) ? val : minv);
00211 maxv = ((val > maxv) ? val : maxv);
00212 }
00213
00214
00215 T step=(maxv-minv)/N;
00216 T w=minv;
00217 for(unsigned int k=0;k<=N;k++,w+=step) d(k)=w;
00218 w=minv+step/2;
00219 for(unsigned int k=0;k<N;k++,w+=step) r(k)=w;
00220
00221 itpp::Vec<T> r_neu(N);
00222 itpp::Vec<unsigned int> r_anz(N);
00223
00224 unsigned int iteration=0;
00225 for(;iteration < max_iter;iteration++)
00226 {
00227 r_neu.clear();
00228 r_anz.clear();
00229
00230 for(int k=0;k<trainingdata.length();k++)
00231 {
00232 T val=trainingdata(k);
00233 T min_dist=std::abs(r(0)-val);
00234 unsigned min_index=0;
00235 for(unsigned int l=1;l<N;l++)
00236 {
00237 T dist=std::abs(T(r(l)-val));
00238 if(dist < min_dist)
00239 {
00240 min_dist=dist;
00241 min_index=l;
00242 }
00243 }
00244 r_neu(min_index)+=val;
00245 r_anz(min_index)++;
00246 }
00247 double change=0.0;
00248 for(unsigned int k=0;k<N;k++)
00249 {
00250 double r_neuw = (r_anz(k) != 0 ? r_neu(k)/r_anz(k) : r(k));
00251 double tmp=r_neuw-r(k);
00252 change+=tmp*tmp;
00253 r(k)=r_neuw;
00254 }
00255
00256
00257
00258 if(change/N <= epsilon) break;
00259 }
00260
00261
00262 for(unsigned int i=1;i<N;i++)
00263 {
00264 T di_neu=(r(i-1)+r(i))/2;
00265
00266 d(i)=di_neu;
00267 }
00268 trained=true;
00269 return (iteration!=max_iter);
00270 }
00271
00272 template<class T> bool MaxLloydQuantizer<T>::train2(const itpp::Vec<T> &trainingdata,double epsilon,unsigned int max_iter)
00273 {
00274 #ifdef DEBUG
00275 if(!initialized) cerr << "MaxLloydQuantizer<T>::train() called before initializing the quantizer!" << endl;
00276 #endif
00277 if(trainingdata.length()==0) return false;
00278 T minv,maxv;
00279 minv=maxv=trainingdata(0);
00280 for(const T* ptr=trainingdata._data();ptr < trainingdata._data()+trainingdata.length();ptr++)
00281 {
00282 register T val = *ptr;
00283 minv = ((val < minv) ? val : minv);
00284 maxv = ((val > maxv) ? val : maxv);
00285 }
00286
00287
00288 T step=(maxv-minv)/N;
00289 T w=minv;
00290 for(unsigned int k=0;k<N;k++,w+=step) d(k)=w;
00291
00292 d(N)=maxv;
00293 w=minv+step/2;
00294 for(unsigned int k=0;k<N;k++,w+=step) r(k)=w;
00295
00296 itpp::Vec<T> r_neu(N);
00297 itpp::Vec<unsigned int> r_anz(N);
00298
00299 unsigned int iteration=0;
00300 for(;iteration < max_iter;iteration++)
00301 {
00302 r_neu.clear();
00303 r_anz.clear();
00304
00305 for(int k=0;k<trainingdata.length();k++)
00306 {
00307 T val=trainingdata(k);
00308 for(unsigned int l=1;l<=N;l++)
00309 {
00310 if(val <= d(l))
00311 {
00312 r_neu(l-1)+=val;
00313 ++r_anz(l-1);
00314 break;
00315 }
00316 }
00317 }
00318 double change=0.0;
00319 for(unsigned int k=0;k<N;k++)
00320 {
00321 r(k) = (r_anz(k) != 0 ? r_neu(k)/r_anz(k) : 0.5*(d(k)+d(k+1)));
00322 }
00323
00324
00325 for (unsigned int i = 1; i < N; i++) {
00326 T di_neu = (r(i - 1) + r(i)) / 2;
00327 change+=std::abs((double) di_neu-d(i));
00328 d(i) = di_neu;
00329 }
00330 cout << "Iteration " << iteration << " var: " << change << endl << " r=" << r << endl;
00331
00332 if(change/N <= epsilon) break;
00333 }
00334
00335 trained=true;
00336 return (iteration!=max_iter);
00337 }
00338
00339
00340 template<class T> unsigned int MaxLloydQuantizer<T>::quantize(T x)
00341 {
00342 #ifdef DEBUG
00343 if(!initialized) cerr << "MaxLloydQuantizer<T>::quant() called before initializing the quantizer!" << endl;
00344 if(!trained) cerr << "MaxLloydQuantizer<T>::quant() called before training the quantizer!" << endl;
00345 #endif
00346 unsigned int k=0;
00347 while(++k != N) if(x < d(k)) break;
00348 return k-1;
00349 }
00350
00351 }
00352 #endif
00353