zeroed channel logic fix to residue backend 2
[platform/upstream/libvorbis.git] / lib / res0.c
1 /********************************************************************
2  *                                                                  *
3  * THIS FILE IS PART OF THE OggVorbis SOFTWARE CODEC SOURCE CODE.   *
4  * USE, DISTRIBUTION AND REPRODUCTION OF THIS LIBRARY SOURCE IS     *
5  * GOVERNED BY A BSD-STYLE SOURCE LICENSE INCLUDED WITH THIS SOURCE *
6  * IN 'COPYING'. PLEASE READ THESE TERMS BEFORE DISTRIBUTING.       *
7  *                                                                  *
8  * THE OggVorbis SOURCE CODE IS (C) COPYRIGHT 1994-2001             *
9  * by the XIPHOPHORUS Company http://www.xiph.org/                  *
10  *                                                                  *
11  ********************************************************************
12
13  function: residue backend 0, 1 and 2 implementation
14  last mod: $Id: res0.c,v 1.32 2001/06/15 23:31:00 xiphmont Exp $
15
16  ********************************************************************/
17
18 /* Slow, slow, slow, simpleminded and did I mention it was slow?  The
19    encode/decode loops are coded for clarity and performance is not
20    yet even a nagging little idea lurking in the shadows.  Oh and BTW,
21    it's slow. */
22
23 #include <stdlib.h>
24 #include <string.h>
25 #include <math.h>
26 #include <stdio.h>
27 #include <ogg/ogg.h>
28 #include "vorbis/codec.h"
29 #include "codec_internal.h"
30 #include "registry.h"
31 #include "codebook.h"
32 #include "misc.h"
33 #include "os.h"
34
35 typedef struct {
36   vorbis_info_residue0 *info;
37   int         map;
38   
39   int         parts;
40   int         stages;
41   codebook   *fullbooks;
42   codebook   *phrasebook;
43   codebook ***partbooks;
44
45   int         partvals;
46   int       **decodemap;
47
48   /*long      resbits[32][32];
49   long      resbitsflat;
50   long      resvals[32];
51   long      phrasebits;
52   long      frames;*/
53
54 } vorbis_look_residue0;
55
56 vorbis_info_residue *res0_copy_info(vorbis_info_residue *vr){
57   vorbis_info_residue0 *info=(vorbis_info_residue0 *)vr;
58   vorbis_info_residue0 *ret=_ogg_malloc(sizeof(vorbis_info_residue0));
59   memcpy(ret,info,sizeof(vorbis_info_residue0));
60   return(ret);
61 }
62
63 void res0_free_info(vorbis_info_residue *i){
64   if(i){
65     memset(i,0,sizeof(vorbis_info_residue0));
66     _ogg_free(i);
67   }
68 }
69
70 void res0_free_look(vorbis_look_residue *i){
71   int j,k;
72   if(i){
73
74     vorbis_look_residue0 *look=(vorbis_look_residue0 *)i;
75     vorbis_info_residue0 *info=look->info;
76
77     /*fprintf(stderr,
78             "%ld frames encoded in %ld phrasebits and %ld residue bits "
79             "(%g/frame) \n",look->frames,look->phrasebits,
80             look->resbitsflat,
81             (look->phrasebits+look->resbitsflat)/(float)look->frames);
82     
83     for(j=0;j<look->parts;j++){
84       long acc=0;
85       fprintf(stderr,"\t[%d] == ",j);
86       for(k=0;k<look->stages;k++)
87         if((info->secondstages[j]>>k)&1){
88           fprintf(stderr,"%ld,",look->resbits[j][k]);
89           acc+=look->resbits[j][k];
90         }
91
92       fprintf(stderr,":: (%ld vals) %1.2fbits/sample\n",look->resvals[j],
93               acc?(float)acc/(look->resvals[j]*info->grouping):0);
94     }
95     fprintf(stderr,"\n");*/
96
97     for(j=0;j<look->parts;j++)
98       if(look->partbooks[j])_ogg_free(look->partbooks[j]);
99     _ogg_free(look->partbooks);
100     for(j=0;j<look->partvals;j++)
101       _ogg_free(look->decodemap[j]);
102     _ogg_free(look->decodemap);
103     memset(i,0,sizeof(vorbis_look_residue0));
104     _ogg_free(i);
105   }
106 }
107
108 static int ilog(unsigned int v){
109   int ret=0;
110   while(v){
111     ret++;
112     v>>=1;
113   }
114   return(ret);
115 }
116
117 static int icount(unsigned int v){
118   int ret=0;
119   while(v){
120     ret+=v&1;
121     v>>=1;
122   }
123   return(ret);
124 }
125
126
127 void res0_pack(vorbis_info_residue *vr,oggpack_buffer *opb){
128   vorbis_info_residue0 *info=(vorbis_info_residue0 *)vr;
129   int j,acc=0;
130   oggpack_write(opb,info->begin,24);
131   oggpack_write(opb,info->end,24);
132
133   oggpack_write(opb,info->grouping-1,24);  /* residue vectors to group and 
134                                              code with a partitioned book */
135   oggpack_write(opb,info->partitions-1,6); /* possible partition choices */
136   oggpack_write(opb,info->groupbook,8);  /* group huffman book */
137
138   /* secondstages is a bitmask; as encoding progresses pass by pass, a
139      bitmask of one indicates this partition class has bits to write
140      this pass */
141   for(j=0;j<info->partitions;j++){
142     if(ilog(info->secondstages[j])>3){
143       /* yes, this is a minor hack due to not thinking ahead */
144       oggpack_write(opb,info->secondstages[j],3); 
145       oggpack_write(opb,1,1);
146       oggpack_write(opb,info->secondstages[j]>>3,5); 
147     }else
148       oggpack_write(opb,info->secondstages[j],4); /* trailing zero */
149     acc+=icount(info->secondstages[j]);
150   }
151   for(j=0;j<acc;j++)
152     oggpack_write(opb,info->booklist[j],8);
153
154 }
155
156 /* vorbis_info is for range checking */
157 vorbis_info_residue *res0_unpack(vorbis_info *vi,oggpack_buffer *opb){
158   int j,acc=0;
159   vorbis_info_residue0 *info=_ogg_calloc(1,sizeof(vorbis_info_residue0));
160   codec_setup_info     *ci=vi->codec_setup;
161
162   info->begin=oggpack_read(opb,24);
163   info->end=oggpack_read(opb,24);
164   info->grouping=oggpack_read(opb,24)+1;
165   info->partitions=oggpack_read(opb,6)+1;
166   info->groupbook=oggpack_read(opb,8);
167
168   for(j=0;j<info->partitions;j++){
169     int cascade=oggpack_read(opb,3);
170     if(oggpack_read(opb,1))
171       cascade|=(oggpack_read(opb,5)<<3);
172     info->secondstages[j]=cascade;
173
174     acc+=icount(cascade);
175   }
176   for(j=0;j<acc;j++)
177     info->booklist[j]=oggpack_read(opb,8);
178
179   if(info->groupbook>=ci->books)goto errout;
180   for(j=0;j<acc;j++)
181     if(info->booklist[j]>=ci->books)goto errout;
182
183   return(info);
184  errout:
185   res0_free_info(info);
186   return(NULL);
187 }
188
189 vorbis_look_residue *res0_look (vorbis_dsp_state *vd,vorbis_info_mode *vm,
190                           vorbis_info_residue *vr){
191   vorbis_info_residue0 *info=(vorbis_info_residue0 *)vr;
192   vorbis_look_residue0 *look=_ogg_calloc(1,sizeof(vorbis_look_residue0));
193   backend_lookup_state *be=vd->backend_state;
194
195   int j,k,acc=0;
196   int dim;
197   int maxstage=0;
198   look->info=info;
199   look->map=vm->mapping;
200
201   look->parts=info->partitions;
202   look->fullbooks=be->fullbooks;
203   look->phrasebook=be->fullbooks+info->groupbook;
204   dim=look->phrasebook->dim;
205
206   look->partbooks=_ogg_calloc(look->parts,sizeof(codebook **));
207
208   for(j=0;j<look->parts;j++){
209     int stages=ilog(info->secondstages[j]);
210     if(stages){
211       if(stages>maxstage)maxstage=stages;
212       look->partbooks[j]=_ogg_calloc(stages,sizeof(codebook *));
213       for(k=0;k<stages;k++)
214         if(info->secondstages[j]&(1<<k))
215           look->partbooks[j][k]=be->fullbooks+info->booklist[acc++];
216     }
217   }
218
219   look->partvals=rint(pow(look->parts,dim));
220   look->stages=maxstage;
221   look->decodemap=_ogg_malloc(look->partvals*sizeof(int *));
222   for(j=0;j<look->partvals;j++){
223     long val=j;
224     long mult=look->partvals/look->parts;
225     look->decodemap[j]=_ogg_malloc(dim*sizeof(int));
226     for(k=0;k<dim;k++){
227       long deco=val/mult;
228       val-=deco*mult;
229       mult/=look->parts;
230       look->decodemap[j][k]=deco;
231     }
232   }
233
234   return(look);
235 }
236
237
238 /* does not guard against invalid settings; eg, a subn of 16 and a
239    subgroup request of 32.  Max subn of 128 */
240 static int _interleaved_testhack(float *vec,int n,vorbis_look_residue0 *look,
241                                  int auxparts,int auxpartnum){
242   vorbis_info_residue0 *info=look->info;
243   int i,j=0;
244   float max,localmax=0.f;
245   float temp[128];
246   float entropy[8];
247
248   /* setup */
249   for(i=0;i<n;i++)temp[i]=fabs(vec[i]);
250
251   /* handle case subgrp==1 outside */
252   for(i=0;i<n;i++)
253     if(temp[i]>localmax)localmax=temp[i];
254   max=localmax;
255
256   for(i=0;i<n;i++)temp[i]=rint(temp[i]);
257   
258   while(1){
259     entropy[j]=localmax;
260     n>>=1;
261     if(!n)break;
262     j++;
263
264     for(i=0;i<n;i++){
265       temp[i]+=temp[i+n];
266     }
267     localmax=0.f;
268     for(i=0;i<n;i++)
269       if(temp[i]>localmax)localmax=temp[i];
270   }
271
272   for(i=0;i<auxparts-1;i++)
273     if(auxpartnum<info->blimit[i] &&
274        entropy[info->subgrp[i]]<=info->entmax[i] &&
275        max<=info->ampmax[i])
276       break;
277
278   return(i);
279 }
280
281 static int _testhack(float *vec,int n,vorbis_look_residue0 *look,
282                      int auxparts,int auxpartnum){
283   vorbis_info_residue0 *info=look->info;
284   int i,j=0;
285   float max,localmax=0.f;
286   float temp[128];
287   float entropy[8];
288
289   /* setup */
290   for(i=0;i<n;i++)temp[i]=fabs(vec[i]);
291
292   /* handle case subgrp==1 outside */
293   for(i=0;i<n;i++)
294     if(temp[i]>localmax)localmax=temp[i];
295   max=localmax;
296
297   for(i=0;i<n;i++)temp[i]=rint(temp[i]);
298   
299   while(n){
300     entropy[j]=localmax;
301     n>>=1;
302     j++;
303     if(!n)break;
304     for(i=0;i<n;i++){
305       temp[i]=temp[i*2]+temp[i*2+1];
306     }
307     localmax=0.f;
308     for(i=0;i<n;i++)
309       if(temp[i]>localmax)localmax=temp[i];
310   }
311
312   for(i=0;i<auxparts-1;i++)
313     if(auxpartnum<info->blimit[i] &&
314        entropy[info->subgrp[i]]<=info->entmax[i] &&
315        max<=info->ampmax[i])
316       break;
317
318   return(i);
319 }
320
321 static int _interleaved_encodepart(oggpack_buffer *opb,float *vec, int n,
322                                    codebook *book,vorbis_look_residue0 *look){
323   int i,bits=0;
324   int dim=book->dim;
325   int step=n/dim;
326 #ifdef TRAIN_RESENT      
327   char buf[80];
328   FILE *f;
329   sprintf(buf,"res0_b%d.vqd",book-look->fullbooks);
330   f=fopen(buf,"a");
331 #endif
332
333   for(i=0;i<step;i++){
334     int entry=vorbis_book_besterror(book,vec+i,step,0);
335
336 #ifdef TRAIN_RESENT      
337     fprintf(f,"%d\n",entry);
338 #endif
339
340     bits+=vorbis_book_encode(book,entry,opb);
341   }
342
343 #ifdef TRAIN_RESENT      
344   fclose(f);
345 #endif
346   return(bits);
347 }
348  
349 static int _encodepart(oggpack_buffer *opb,float *vec, int n,
350                        codebook *book,vorbis_look_residue0 *look){
351   int i,bits=0;
352   int dim=book->dim;
353   int step=n/dim;
354 #ifdef TRAIN_RESENT      
355   char buf[80];
356   FILE *f;
357   sprintf(buf,"res0_b%d.vqd",book-look->fullbooks);
358   f=fopen(buf,"a");
359 #endif
360
361   for(i=0;i<step;i++){
362     int entry=vorbis_book_besterror(book,vec+i*dim,1,0);
363
364 #ifdef TRAIN_RESENT      
365     fprintf(f,"%d\n",entry);
366 #endif
367
368     bits+=vorbis_book_encode(book,entry,opb);
369   }
370
371 #ifdef TRAIN_RESENT      
372   fclose(f);
373 #endif
374   return(bits);
375 }
376
377 static int _01forward(vorbis_block *vb,vorbis_look_residue *vl,
378                       float **in,int ch,
379                       int (*classify)(float *,int,vorbis_look_residue0 *,
380                                       int,int),
381                       int (*encode)(oggpack_buffer *,float *,int,
382                                     codebook *,vorbis_look_residue0 *)){
383   long i,j,k,l,s;
384   vorbis_look_residue0 *look=(vorbis_look_residue0 *)vl;
385   vorbis_info_residue0 *info=look->info;
386
387   /* move all this setup out later */
388   int samples_per_partition=info->grouping;
389   int possible_partitions=info->partitions;
390   int partitions_per_word=look->phrasebook->dim;
391   int n=info->end-info->begin;
392
393   int partvals=n/samples_per_partition;
394   int partwords=(partvals+partitions_per_word-1)/partitions_per_word;
395   long **partword=_vorbis_block_alloc(vb,ch*sizeof(long *));
396
397 #ifdef TRAIN_RES
398   FILE *of;
399   char buffer[80];
400   int m;
401   
402   for(i=0;i<ch;i++){
403     sprintf(buffer,"residue_%d.vqd",vb->mode);
404     of=fopen(buffer,"a");
405     for(m=0;m<info->end;m++)
406       fprintf(of,"%.2f, ",in[i][m]);
407     fprintf(of,"\n");
408     fclose(of);
409   }
410  
411 #endif      
412
413   partvals=partwords*partitions_per_word;
414
415   /* we find the patition type for each partition of each
416      channel.  We'll go back and do the interleaved encoding in a
417      bit.  For now, clarity */
418  
419   for(i=0;i<ch;i++){
420     partword[i]=_vorbis_block_alloc(vb,n/samples_per_partition*sizeof(long));
421     memset(partword[i],0,n/samples_per_partition*sizeof(long));
422   }
423
424   for(i=info->begin,l=0;i<info->end;i+=samples_per_partition,l++){
425     for(j=0;j<ch;j++)
426       /* do the partition decision based on the 'entropy'
427          int the block */
428       partword[j][l]=
429         classify(in[j]+i,samples_per_partition,look,possible_partitions,l);
430   
431   }
432
433   /* we code the partition words for each channel, then the residual
434      words for a partition per channel until we've written all the
435      residual words for that partition word.  Then write the next
436      partition channel words... */
437
438   /*look->frames++;*/
439   for(s=0;s<look->stages;s++){
440     for(i=info->begin,l=0;i<info->end;){
441       
442       /* first we encode a partition codeword for each channel */
443       if(s==0){
444         for(j=0;j<ch;j++){
445           long val=partword[j][l];
446           long ret;
447           for(k=1;k<partitions_per_word;k++)
448             val= val*possible_partitions+partword[j][l+k];
449           ret=vorbis_book_encode(look->phrasebook,val,&vb->opb);
450           /*look->phrasebits+=ret;*/
451         }
452       }
453       
454       /* now we encode interleaved residual values for the partitions */
455       for(k=0;k<partitions_per_word;k++,l++,i+=samples_per_partition){
456         
457         for(j=0;j<ch;j++){
458           /*if(s==0)look->resvals[partword[j][l]]++;*/
459           if(info->secondstages[partword[j][l]]&(1<<s)){
460             codebook *statebook=look->partbooks[partword[j][l]][s];
461             if(statebook){
462               int ret=encode(&vb->opb,in[j]+i,samples_per_partition,
463                              statebook,look);
464               /*look->resbits[partword[j][l]][s]+=ret;
465                 look->resbitsflat+=ret;*/
466
467             }
468           }
469         }
470       }
471     }
472   }
473   
474   return(0);
475 }
476
477 /* a truncated packet here just means 'stop working'; it's not an error */
478 static int _01inverse(vorbis_block *vb,vorbis_look_residue *vl,
479                       float **in,int ch,
480                       long (*decodepart)(codebook *, float *, 
481                                          oggpack_buffer *,int)){
482
483   long i,j,k,l,s;
484   vorbis_look_residue0 *look=(vorbis_look_residue0 *)vl;
485   vorbis_info_residue0 *info=look->info;
486
487   /* move all this setup out later */
488   int samples_per_partition=info->grouping;
489   int partitions_per_word=look->phrasebook->dim;
490   int n=info->end-info->begin;
491
492   int partvals=n/samples_per_partition;
493   int partwords=(partvals+partitions_per_word-1)/partitions_per_word;
494   int ***partword=alloca(ch*sizeof(int **));
495   partvals=partwords*partitions_per_word;
496
497   for(j=0;j<ch;j++)
498     partword[j]=_vorbis_block_alloc(vb,partwords*sizeof(int *));
499
500   for(s=0;s<look->stages;s++){
501     for(i=info->begin,l=0;i<info->end;l++){
502
503       if(s==0){
504         /* fetch the partition word for each channel */
505         for(j=0;j<ch;j++){
506           int temp=vorbis_book_decode(look->phrasebook,&vb->opb);
507           if(temp==-1)goto eopbreak;
508           partword[j][l]=look->decodemap[temp];
509           if(partword[j][l]==NULL)goto errout;
510         }
511       }
512       
513       /* now we decode residual values for the partitions */
514       for(k=0;k<partitions_per_word;k++,i+=samples_per_partition)
515         for(j=0;j<ch;j++){
516           if(info->secondstages[partword[j][l][k]]&(1<<s)){
517             codebook *stagebook=look->partbooks[partword[j][l][k]][s];
518             if(stagebook){
519               if(decodepart(stagebook,in[j]+i,&vb->opb,
520                             samples_per_partition)==-1)goto eopbreak;
521             }
522           }
523         }
524     } 
525   }
526   
527  errout:
528  eopbreak:
529   return(0);
530 }
531
532 /* residue 0 and 1 are just slight variants of one another. 0 is
533    interleaved, 1 is not */
534 int res0_forward(vorbis_block *vb,vorbis_look_residue *vl,
535             float **in,int *nonzero,int ch){
536   /* we encode only the nonzero parts of a bundle */
537   int i,used=0;
538   for(i=0;i<ch;i++)
539     if(nonzero[i])
540       in[used++]=in[i];
541   if(used)
542     return(_01forward(vb,vl,in,used,_interleaved_testhack,_interleaved_encodepart));
543   else
544     return(0);
545 }
546
547 int res0_inverse(vorbis_block *vb,vorbis_look_residue *vl,
548                  float **in,int *nonzero,int ch){
549   int i,used=0;
550   for(i=0;i<ch;i++)
551     if(nonzero[i])
552       in[used++]=in[i];
553   if(used)
554     return(_01inverse(vb,vl,in,used,vorbis_book_decodevs_add));
555   else
556     return(0);
557 }
558
559 int res1_forward(vorbis_block *vb,vorbis_look_residue *vl,
560                  float **in,int *nonzero,int ch){
561   int i,used=0;
562   for(i=0;i<ch;i++)
563     if(nonzero[i])
564       in[used++]=in[i];
565   if(used)
566     return(_01forward(vb,vl,in,used,_testhack,_encodepart));
567   else
568     return(0);
569 }
570
571 int res1_inverse(vorbis_block *vb,vorbis_look_residue *vl,
572                  float **in,int *nonzero,int ch){
573   int i,used=0;
574   for(i=0;i<ch;i++)
575     if(nonzero[i])
576       in[used++]=in[i];
577   if(used)
578     return(_01inverse(vb,vl,in,used,vorbis_book_decodev_add));
579   else
580     return(0);
581 }
582
583 /* res2 is slightly more different; all the channels are interleaved
584    into a single vector and encoded. */
585 int res2_forward(vorbis_block *vb,vorbis_look_residue *vl,
586             float **in,int *nonzero,int ch){
587   long i,j,k,n=vb->pcmend/2,used=0;
588
589   /* don't duplicate the code; use a working vector hack for now and
590      reshape ourselves into a single channel res1 */
591   float *work=_vorbis_block_alloc(vb,ch*n*sizeof(float));
592   for(i=0;i<ch;i++){
593     float *pcm=vb->pcm[i];
594     if(nonzero[i])used++;
595     for(j=0,k=i;j<n;j++,k+=ch)
596       work[k]=pcm[j];
597   }
598   
599   if(used)
600     return(_01forward(vb,vl,&work,1,_testhack,_encodepart));
601   else
602     return(0);
603 }
604
605 /* duplicate code here as speed is somewhat more important */
606 int res2_inverse(vorbis_block *vb,vorbis_look_residue *vl,
607                  float **in,int *nonzero,int ch){
608   long i,k,l,s;
609   vorbis_look_residue0 *look=(vorbis_look_residue0 *)vl;
610   vorbis_info_residue0 *info=look->info;
611
612   /* move all this setup out later */
613   int samples_per_partition=info->grouping;
614   int partitions_per_word=look->phrasebook->dim;
615   int n=info->end-info->begin;
616
617   int partvals=n/samples_per_partition;
618   int partwords=(partvals+partitions_per_word-1)/partitions_per_word;
619   int **partword=_vorbis_block_alloc(vb,partwords*sizeof(int *));
620   int used;
621   partvals=partwords*partitions_per_word;
622
623   for(i=0;i<ch;i++)if(nonzero[i])break;
624   if(i==ch)return(0); /* no nonzero vectors */
625
626   for(s=0;s<look->stages;s++){
627     for(i=info->begin,l=0;i<info->end;l++){
628
629       if(s==0){
630         /* fetch the partition word */
631         int temp=vorbis_book_decode(look->phrasebook,&vb->opb);
632         if(temp==-1)goto eopbreak;
633         partword[l]=look->decodemap[temp];
634         if(partword[l]==NULL)goto errout;
635       }
636
637       /* now we decode residual values for the partitions */
638       for(k=0;k<partitions_per_word;k++,i+=samples_per_partition)
639         if(info->secondstages[partword[l][k]]&(1<<s)){
640           codebook *stagebook=look->partbooks[partword[l][k]][s];
641
642           if(stagebook){
643             if(vorbis_book_decodevv_add(stagebook,in,i,ch,
644                                         &vb->opb,samples_per_partition)==-1)
645               goto eopbreak;
646           }
647         }
648     } 
649   }
650   
651  errout:
652  eopbreak:
653   return(0);
654 }
655
656
657 vorbis_func_residue residue0_exportbundle={
658   &res0_pack,
659   &res0_unpack,
660   &res0_look,
661   &res0_copy_info,
662   &res0_free_info,
663   &res0_free_look,
664   &res0_forward,
665   &res0_inverse
666 };
667
668 vorbis_func_residue residue1_exportbundle={
669   &res0_pack,
670   &res0_unpack,
671   &res0_look,
672   &res0_copy_info,
673   &res0_free_info,
674   &res0_free_look,
675   &res1_forward,
676   &res1_inverse
677 };
678
679 vorbis_func_residue residue2_exportbundle={
680   &res0_pack,
681   &res0_unpack,
682   &res0_look,
683   &res0_copy_info,
684   &res0_free_info,
685   &res0_free_look,
686   &res2_forward,
687   &res2_inverse
688 };