envelope.c 6.0 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220
  1. /********************************************************************
  2. * *
  3. * THIS FILE IS PART OF THE OggVorbis SOFTWARE CODEC SOURCE CODE. *
  4. * USE, DISTRIBUTION AND REPRODUCTION OF THIS SOURCE IS GOVERNED BY *
  5. * THE GNU LESSER/LIBRARY PUBLIC LICENSE, WHICH IS INCLUDED WITH *
  6. * THIS SOURCE. PLEASE READ THESE TERMS BEFORE DISTRIBUTING. *
  7. * *
  8. * THE OggVorbis SOURCE CODE IS (C) COPYRIGHT 1994-2000 *
  9. * by Monty <monty@xiph.org> and the XIPHOPHORUS Company *
  10. * http://www.xiph.org/ *
  11. * *
  12. ********************************************************************
  13. function: PCM data envelope analysis and manipulation
  14. last mod: $Id: envelope.c,v 1.23.2.4 2000/11/04 10:24:15 xiphmont Exp $
  15. Preecho calculation.
  16. ********************************************************************/
  17. #include <stdlib.h>
  18. #include <string.h>
  19. #include <stdio.h>
  20. #include <math.h>
  21. #include <ogg/ogg.h>
  22. #include "vorbis/codec.h"
  23. #include "codec_internal.h"
  24. #include "os.h"
  25. #include "scales.h"
  26. #include "envelope.h"
  27. #include "misc.h"
  28. /* We use a Chebyshev bandbass for the preecho trigger bandpass; it's
  29. close enough for sample rates 32000-48000 Hz (corner frequencies at
  30. 6k/14k assuming sample rate of 44.1kHz) */
  31. /* Digital filter designed by mkfilter/mkshape/gencode A.J. Fisher
  32. Command line: /www/usr/fisher/helpers/mkfilter -Ch \
  33. -6.0000000000e+00 -Bp -o 5 -a 1.3605442177e-01 3.1746031746e-01 -l */
  34. #if 0
  35. static int cheb_bandpass_stages=10;
  36. static float cheb_bandpass_gain=5.589612458e+01;
  37. static float cheb_bandpass_B[]={-1.,0.,5.,0.,-10.,0.,10.,0.,-5.,0.,1};
  38. static float cheb_bandpass_A[]={
  39. -0.1917409386,
  40. 0.0078657069,
  41. -0.7126903444,
  42. 0.0266343467,
  43. -1.4047174730,
  44. 0.0466964232,
  45. -1.9032773429,
  46. 0.0451493360,
  47. -1.4471447397,
  48. 0.0303413711};
  49. #endif
  50. static int cheb_highpass_stages=10;
  51. static float cheb_highpass_gain= 5.291963434e+01;
  52. /* z^-stage, z^-stage+1... */
  53. static float cheb_highpass_B[]={1,-10,45,-120,210,-252,210,-120,45,-10,1};
  54. static float cheb_highpass_A[]={
  55. -0.1247628029,
  56. 0.1334086523,
  57. -0.3997715614,
  58. 0.3213011089,
  59. -1.1131924119,
  60. 1.7692446626,
  61. -3.6241199038,
  62. 4.1950871291,
  63. -4.2771757867,
  64. 2.3920318913};
  65. void _ve_envelope_init(envelope_lookup *e,vorbis_info *vi){
  66. codec_setup_info *ci=vi->codec_setup;
  67. int ch=vi->channels;
  68. int window=ci->envelopesa;
  69. int i;
  70. e->winlength=window;
  71. e->minenergy=fromdB(ci->preecho_minenergy);
  72. e->iir=_ogg_calloc(ch,sizeof(IIR_state));
  73. e->filtered=_ogg_calloc(ch,sizeof(float *));
  74. e->ch=ch;
  75. e->storage=128;
  76. for(i=0;i<ch;i++){
  77. IIR_init(e->iir+i,cheb_highpass_stages,cheb_highpass_gain,
  78. cheb_highpass_A,cheb_highpass_B);
  79. e->filtered[i]=_ogg_calloc(e->storage,sizeof(float));
  80. }
  81. drft_init(&e->drft,window);
  82. e->window=_ogg_malloc(e->winlength*sizeof(float));
  83. /* We just use a straight sin(x) window for this */
  84. for(i=0;i<e->winlength;i++)
  85. e->window[i]=sin((i+.5)/e->winlength*M_PI);
  86. }
  87. void _ve_envelope_clear(envelope_lookup *e){
  88. int i;
  89. for(i=0;i<e->ch;i++){
  90. IIR_clear((e->iir+i));
  91. free(e->filtered[i]);
  92. }
  93. drft_clear(&e->drft);
  94. free(e->window);
  95. free(e->filtered);
  96. memset(e,0,sizeof(envelope_lookup));
  97. }
  98. static float _ve_deltai(envelope_lookup *ve,IIR_state *iir,
  99. float *pre,float *post){
  100. long n2=ve->winlength*2;
  101. long n=ve->winlength;
  102. float *workA=alloca(sizeof(float)*n2),A=0.;
  103. float *workB=alloca(sizeof(float)*n2),B=0.;
  104. long i;
  105. /*_analysis_output("A",granulepos,pre,n,0,0);
  106. _analysis_output("B",granulepos,post,n,0,0);*/
  107. for(i=0;i<n;i++){
  108. workA[i]=pre[i]*ve->window[i];
  109. workB[i]=post[i]*ve->window[i];
  110. }
  111. /*_analysis_output("Awin",granulepos,workA,n,0,0);
  112. _analysis_output("Bwin",granulepos,workB,n,0,0);*/
  113. drft_forward(&ve->drft,workA);
  114. drft_forward(&ve->drft,workB);
  115. /* we want to have a 'minimum bar' for energy, else we're just
  116. basing blocks on quantization noise that outweighs the signal
  117. itself (for low power signals) */
  118. {
  119. float min=ve->minenergy;
  120. for(i=0;i<n;i++){
  121. if(fabs(workA[i])<min)workA[i]=min;
  122. if(fabs(workB[i])<min)workB[i]=min;
  123. }
  124. }
  125. /*_analysis_output("Afft",granulepos,workA,n,0,0);
  126. _analysis_output("Bfft",granulepos,workB,n,0,0);*/
  127. for(i=0;i<n;i++){
  128. A+=workA[i]*workA[i];
  129. B+=workB[i]*workB[i];
  130. }
  131. A=todB(A);
  132. B=todB(B);
  133. return(B-A);
  134. }
  135. long _ve_envelope_search(vorbis_dsp_state *v,long searchpoint){
  136. vorbis_info *vi=v->vi;
  137. codec_setup_info *ci=vi->codec_setup;
  138. envelope_lookup *ve=((backend_lookup_state *)(v->backend_state))->ve;
  139. long i,j;
  140. /* make sure we have enough storage to match the PCM */
  141. if(v->pcm_storage>ve->storage){
  142. ve->storage=v->pcm_storage;
  143. for(i=0;i<ve->ch;i++)
  144. ve->filtered[i]=_ogg_realloc(ve->filtered[i],ve->storage*sizeof(float));
  145. }
  146. /* catch up the highpass to match the pcm */
  147. for(i=0;i<ve->ch;i++){
  148. float *filtered=ve->filtered[i];
  149. float *pcm=v->pcm[i];
  150. IIR_state *iir=ve->iir+i;
  151. IIR_clamp(iir,9e-15);
  152. for(j=ve->current;j<v->pcm_current;j++)
  153. filtered[j]=IIR_filter(iir,pcm[j]);
  154. }
  155. ve->current=v->pcm_current;
  156. /* Now search through our cached highpass data for breaking points */
  157. /* starting point */
  158. if(v->W)
  159. j=v->centerW+ci->blocksizes[1]/4-ci->blocksizes[0]/4;
  160. else
  161. j=v->centerW;
  162. while(j+ve->winlength<=v->pcm_current){
  163. for(i=0;i<ve->ch;i++){
  164. float *filtered=ve->filtered[i]+j;
  165. IIR_state *iir=ve->iir+i;
  166. float m=_ve_deltai(ve,iir,filtered-ve->winlength,filtered);
  167. if(m>ci->preecho_thresh){
  168. /*granulepos++;*/
  169. return(0);
  170. }
  171. /*granulepos++;*/
  172. }
  173. j+=ci->blocksizes[0]/2;
  174. if(j>=searchpoint)return(1);
  175. }
  176. return(-1);
  177. }
  178. void _ve_envelope_shift(envelope_lookup *e,long shift){
  179. int i;
  180. for(i=0;i<e->ch;i++)
  181. memmove(e->filtered[i],e->filtered[i]+shift,(e->current-shift)*
  182. sizeof(float));
  183. e->current-=shift;
  184. }