floor0.c 12 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448
  1. /************************************************************************
  2. * Copyright (C) 2002-2009, Xiph.org Foundation
  3. * Copyright (C) 2010, Robin Watts for Pinknoise Productions Ltd
  4. * All rights reserved.
  5. *
  6. * Redistribution and use in source and binary forms, with or without
  7. * modification, are permitted provided that the following conditions
  8. * are met:
  9. *
  10. * * Redistributions of source code must retain the above copyright
  11. * notice, this list of conditions and the following disclaimer.
  12. * * Redistributions in binary form must reproduce the above
  13. * copyright notice, this list of conditions and the following disclaimer
  14. * in the documentation and/or other materials provided with the
  15. * distribution.
  16. * * Neither the names of the Xiph.org Foundation nor Pinknoise
  17. * Productions Ltd nor the names of its contributors may be used to
  18. * endorse or promote products derived from this software without
  19. * specific prior written permission.
  20. *
  21. * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
  22. * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
  23. * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
  24. * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
  25. * OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
  26. * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
  27. * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
  28. * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
  29. * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
  30. * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
  31. * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
  32. ************************************************************************
  33. function: floor backend 0 implementation
  34. ************************************************************************/
  35. #include <stdlib.h>
  36. #include <string.h>
  37. #include <math.h>
  38. #include "ogg.h"
  39. #include "ivorbiscodec.h"
  40. #include "codec_internal.h"
  41. #include "codebook.h"
  42. #include "misc.h"
  43. #include "os.h"
  44. #define LSP_FRACBITS 14
  45. extern const ogg_int32_t FLOOR_fromdB_LOOKUP[];
  46. /*************** LSP decode ********************/
  47. #include "lsp_lookup.h"
  48. /* interpolated 1./sqrt(p) where .5 <= a < 1. (.100000... to .111111...) in
  49. 16.16 format
  50. returns in m.8 format */
  51. static long ADJUST_SQRT2[2]={8192,5792};
  52. static inline ogg_int32_t vorbis_invsqlook_i(long a,long e){
  53. long i=(a&0x7fff)>>(INVSQ_LOOKUP_I_SHIFT-1);
  54. long d=a&INVSQ_LOOKUP_I_MASK; /* 0.10 */
  55. long val=INVSQ_LOOKUP_I[i]- /* 1.16 */
  56. ((INVSQ_LOOKUP_IDel[i]*d)>>INVSQ_LOOKUP_I_SHIFT); /* result 1.16 */
  57. val*=ADJUST_SQRT2[e&1];
  58. e=(e>>1)+21;
  59. return(val>>e);
  60. }
  61. /* interpolated lookup based fromdB function, domain -140dB to 0dB only */
  62. /* a is in n.12 format */
  63. #ifdef _LOW_ACCURACY_
  64. static inline ogg_int32_t vorbis_fromdBlook_i(long a){
  65. if(a>0) return 0x7fffffff;
  66. if(a<(int)(((unsigned)-140)<<12)) return 0;
  67. return FLOOR_fromdB_LOOKUP[((a+140)*467)>>20]<<9;
  68. }
  69. #else
  70. static inline ogg_int32_t vorbis_fromdBlook_i(long a){
  71. if(a>0) return 0x7fffffff;
  72. if(a<(int)(((unsigned)-140)<<12)) return 0;
  73. return FLOOR_fromdB_LOOKUP[((a+(140<<12))*467)>>20];
  74. }
  75. #endif
  76. /* interpolated lookup based cos function, domain 0 to PI only */
  77. /* a is in 0.16 format, where 0==0, 2^^16-1==PI, return 0.14 */
  78. static inline ogg_int32_t vorbis_coslook_i(long a){
  79. int i=a>>COS_LOOKUP_I_SHIFT;
  80. int d=a&COS_LOOKUP_I_MASK;
  81. return COS_LOOKUP_I[i]- ((d*(COS_LOOKUP_I[i]-COS_LOOKUP_I[i+1]))>>
  82. COS_LOOKUP_I_SHIFT);
  83. }
  84. /* interpolated half-wave lookup based cos function */
  85. /* a is in 0.16 format, where 0==0, 2^^16==PI, return .LSP_FRACBITS */
  86. static inline ogg_int32_t vorbis_coslook2_i(long a){
  87. int i=a>>COS_LOOKUP_I_SHIFT;
  88. int d=a&COS_LOOKUP_I_MASK;
  89. return ((COS_LOOKUP_I[i]<<COS_LOOKUP_I_SHIFT)-
  90. d*(COS_LOOKUP_I[i]-COS_LOOKUP_I[i+1]))>>
  91. (COS_LOOKUP_I_SHIFT-LSP_FRACBITS+14);
  92. }
  93. static const ogg_uint16_t barklook[54]={
  94. 0,51,102,154, 206,258,311,365,
  95. 420,477,535,594, 656,719,785,854,
  96. 926,1002,1082,1166, 1256,1352,1454,1564,
  97. 1683,1812,1953,2107, 2276,2463,2670,2900,
  98. 3155,3440,3756,4106, 4493,4919,5387,5901,
  99. 6466,7094,7798,8599, 9528,10623,11935,13524,
  100. 15453,17775,20517,23667, 27183,31004
  101. };
  102. /* used in init only; interpolate the long way */
  103. static inline ogg_int32_t toBARK(int n){
  104. int i;
  105. for(i=0;i<54;i++)
  106. if(n>=barklook[i] && n<barklook[i+1])break;
  107. if(i==54){
  108. return 54<<14;
  109. }else{
  110. return (i<<14)+(((n-barklook[i])*
  111. ((1UL<<31)/(barklook[i+1]-barklook[i])))>>17);
  112. }
  113. }
  114. static const unsigned char MLOOP_1[64]={
  115. 0,10,11,11, 12,12,12,12, 13,13,13,13, 13,13,13,13,
  116. 14,14,14,14, 14,14,14,14, 14,14,14,14, 14,14,14,14,
  117. 15,15,15,15, 15,15,15,15, 15,15,15,15, 15,15,15,15,
  118. 15,15,15,15, 15,15,15,15, 15,15,15,15, 15,15,15,15,
  119. };
  120. static const unsigned char MLOOP_2[64]={
  121. 0,4,5,5, 6,6,6,6, 7,7,7,7, 7,7,7,7,
  122. 8,8,8,8, 8,8,8,8, 8,8,8,8, 8,8,8,8,
  123. 9,9,9,9, 9,9,9,9, 9,9,9,9, 9,9,9,9,
  124. 9,9,9,9, 9,9,9,9, 9,9,9,9, 9,9,9,9,
  125. };
  126. static const unsigned char MLOOP_3[8]={0,1,2,2,3,3,3,3};
  127. void vorbis_lsp_to_curve(ogg_int32_t *curve,int n,int ln,
  128. ogg_int32_t *lsp,int m,
  129. ogg_int32_t amp,
  130. ogg_int32_t ampoffset,
  131. ogg_int32_t nyq){
  132. /* 0 <= m < 256 */
  133. /* set up for using all int later */
  134. int i;
  135. int ampoffseti=ampoffset*4096;
  136. int ampi=amp;
  137. ogg_int32_t *ilsp=(ogg_int32_t *)alloca(m*sizeof(*ilsp));
  138. ogg_uint32_t inyq= (1UL<<31) / toBARK(nyq);
  139. ogg_uint32_t imap= (1UL<<31) / ln;
  140. ogg_uint32_t tBnyq1 = toBARK(nyq)<<1;
  141. /* Besenham for frequency scale to avoid a division */
  142. int f=0;
  143. int fdx=n;
  144. int fbase=nyq/fdx;
  145. int ferr=0;
  146. int fdy=nyq-fbase*fdx;
  147. int map=0;
  148. #ifdef _LOW_ACCURACY_
  149. ogg_uint32_t nextbark=((tBnyq1<<11)/ln)>>12;
  150. #else
  151. ogg_uint32_t nextbark=MULT31(imap>>1,tBnyq1);
  152. #endif
  153. int nextf=barklook[nextbark>>14]+(((nextbark&0x3fff)*
  154. (barklook[(nextbark>>14)+1]-barklook[nextbark>>14]))>>14);
  155. /* lsp is in 8.24, range 0 to PI; coslook wants it in .16 0 to 1*/
  156. for(i=0;i<m;i++){
  157. #ifndef _LOW_ACCURACY_
  158. ogg_int32_t val=MULT32(lsp[i],0x517cc2);
  159. #else
  160. ogg_int32_t val=((lsp[i]>>10)*0x517d)>>14;
  161. #endif
  162. /* safeguard against a malicious stream */
  163. if(val<0 || (val>>COS_LOOKUP_I_SHIFT)>=COS_LOOKUP_I_SZ){
  164. memset(curve,0,sizeof(*curve)*n);
  165. return;
  166. }
  167. ilsp[i]=vorbis_coslook_i(val);
  168. }
  169. i=0;
  170. while(i<n){
  171. int j;
  172. ogg_uint32_t pi=46341; /* 2**-.5 in 0.16 */
  173. ogg_uint32_t qi=46341;
  174. ogg_int32_t qexp=0,shift;
  175. ogg_int32_t wi;
  176. wi=vorbis_coslook2_i((map*imap)>>15);
  177. #ifdef _V_LSP_MATH_ASM
  178. lsp_loop_asm(&qi,&pi,&qexp,ilsp,wi,m);
  179. pi=((pi*pi)>>16);
  180. qi=((qi*qi)>>16);
  181. if(m&1){
  182. qexp= qexp*2-28*((m+1)>>1)+m;
  183. pi*=(1<<14)-((wi*wi)>>14);
  184. qi+=pi>>14;
  185. }else{
  186. qexp= qexp*2-13*m;
  187. pi*=(1<<14)-wi;
  188. qi*=(1<<14)+wi;
  189. qi=(qi+pi)>>14;
  190. }
  191. if(qi&0xffff0000){ /* checks for 1.xxxxxxxxxxxxxxxx */
  192. qi>>=1; qexp++;
  193. }else
  194. lsp_norm_asm(&qi,&qexp);
  195. #else
  196. qi*=labs(ilsp[0]-wi);
  197. pi*=labs(ilsp[1]-wi);
  198. for(j=3;j<m;j+=2){
  199. if(!(shift=MLOOP_1[(pi|qi)>>25]))
  200. if(!(shift=MLOOP_2[(pi|qi)>>19]))
  201. shift=MLOOP_3[(pi|qi)>>16];
  202. qi=(qi>>shift)*labs(ilsp[j-1]-wi);
  203. pi=(pi>>shift)*labs(ilsp[j]-wi);
  204. qexp+=shift;
  205. }
  206. if(!(shift=MLOOP_1[(pi|qi)>>25]))
  207. if(!(shift=MLOOP_2[(pi|qi)>>19]))
  208. shift=MLOOP_3[(pi|qi)>>16];
  209. /* pi,qi normalized collectively, both tracked using qexp */
  210. if(m&1){
  211. /* odd order filter; slightly assymetric */
  212. /* the last coefficient */
  213. qi=(qi>>shift)*labs(ilsp[j-1]-wi);
  214. pi=(pi>>shift)<<14;
  215. qexp+=shift;
  216. if(!(shift=MLOOP_1[(pi|qi)>>25]))
  217. if(!(shift=MLOOP_2[(pi|qi)>>19]))
  218. shift=MLOOP_3[(pi|qi)>>16];
  219. pi>>=shift;
  220. qi>>=shift;
  221. qexp+=shift-14*((m+1)>>1);
  222. pi=((pi*pi)>>16);
  223. qi=((qi*qi)>>16);
  224. qexp=qexp*2+m;
  225. pi*=(1<<14)-((wi*wi)>>14);
  226. qi+=pi>>14;
  227. }else{
  228. /* even order filter; still symmetric */
  229. /* p*=p(1-w), q*=q(1+w), let normalization drift because it isn't
  230. worth tracking step by step */
  231. pi>>=shift;
  232. qi>>=shift;
  233. qexp+=shift-7*m;
  234. pi=((pi*pi)>>16);
  235. qi=((qi*qi)>>16);
  236. qexp=qexp*2+m;
  237. pi*=(1<<14)-wi;
  238. qi*=(1<<14)+wi;
  239. qi=(qi+pi)>>14;
  240. }
  241. /* we've let the normalization drift because it wasn't important;
  242. however, for the lookup, things must be normalized again. We
  243. need at most one right shift or a number of left shifts */
  244. if(qi&0xffff0000){ /* checks for 1.xxxxxxxxxxxxxxxx */
  245. qi>>=1; qexp++;
  246. }else
  247. while(qi && !(qi&0x8000)){ /* checks for 0.0xxxxxxxxxxxxxxx or less*/
  248. qi<<=1; qexp--;
  249. }
  250. #endif
  251. amp=vorbis_fromdBlook_i(ampi* /* n.4 */
  252. vorbis_invsqlook_i(qi,qexp)-
  253. /* m.8, m+n<=8 */
  254. ampoffseti); /* 8.12[0] */
  255. #ifdef _LOW_ACCURACY_
  256. amp>>=9;
  257. #endif
  258. curve[i]= MULT31_SHIFT15(curve[i],amp);
  259. while(++i<n){
  260. /* line plot to get new f */
  261. ferr+=fdy;
  262. if(ferr>=fdx){
  263. ferr-=fdx;
  264. f++;
  265. }
  266. f+=fbase;
  267. if(f>=nextf)break;
  268. curve[i]= MULT31_SHIFT15(curve[i],amp);
  269. }
  270. while(1){
  271. map++;
  272. if(map+1<ln){
  273. #ifdef _LOW_ACCURACY_
  274. nextbark=((tBnyq1<<11)/ln*(map+1))>>12;
  275. #else
  276. nextbark=MULT31((map+1)*(imap>>1),tBnyq1);
  277. #endif
  278. nextf=barklook[nextbark>>14]+
  279. (((nextbark&0x3fff)*
  280. (barklook[(nextbark>>14)+1]-barklook[nextbark>>14]))>>14);
  281. if(f<=nextf)break;
  282. }else{
  283. nextf=9999999;
  284. break;
  285. }
  286. }
  287. if(map>=ln){
  288. map=ln-1; /* guard against the approximation */
  289. nextf=9999999;
  290. }
  291. }
  292. }
  293. /*************** vorbis decode glue ************/
  294. void floor0_free_info(vorbis_info_floor *i){
  295. vorbis_info_floor0 *info=(vorbis_info_floor0 *)i;
  296. if(info)_ogg_free(info);
  297. }
  298. vorbis_info_floor *floor0_info_unpack (vorbis_info *vi,oggpack_buffer *opb){
  299. codec_setup_info *ci=(codec_setup_info *)vi->codec_setup;
  300. int j;
  301. vorbis_info_floor0 *info=(vorbis_info_floor0 *)_ogg_malloc(sizeof(*info));
  302. info->order=oggpack_read(opb,8);
  303. info->rate=oggpack_read(opb,16);
  304. info->barkmap=oggpack_read(opb,16);
  305. info->ampbits=oggpack_read(opb,6);
  306. info->ampdB=oggpack_read(opb,8);
  307. info->numbooks=oggpack_read(opb,4)+1;
  308. if(info->order<1)goto err_out;
  309. if(info->rate<1)goto err_out;
  310. if(info->barkmap<1)goto err_out;
  311. for(j=0;j<info->numbooks;j++){
  312. info->books[j]=(char)oggpack_read(opb,8);
  313. if(info->books[j]>=ci->books)goto err_out;
  314. }
  315. if(oggpack_eop(opb))goto err_out;
  316. return(info);
  317. err_out:
  318. floor0_free_info(info);
  319. return(NULL);
  320. }
  321. int floor0_memosize(vorbis_info_floor *i){
  322. vorbis_info_floor0 *info=(vorbis_info_floor0 *)i;
  323. return info->order+1;
  324. }
  325. ogg_int32_t *floor0_inverse1(vorbis_dsp_state *vd,vorbis_info_floor *i,
  326. ogg_int32_t *lsp){
  327. vorbis_info_floor0 *info=(vorbis_info_floor0 *)i;
  328. int j,k;
  329. int ampraw=oggpack_read(&vd->opb,info->ampbits);
  330. if(ampraw>0){ /* also handles the -1 out of data case */
  331. long maxval=(1<<info->ampbits)-1;
  332. int amp=((ampraw*info->ampdB)<<4)/maxval;
  333. int booknum=oggpack_read(&vd->opb,_ilog(info->numbooks));
  334. if(booknum!=-1 && booknum<info->numbooks){ /* be paranoid */
  335. codec_setup_info *ci=(codec_setup_info *)vd->vi->codec_setup;
  336. codebook *b=ci->book_param+info->books[booknum];
  337. ogg_int32_t last=0;
  338. for(j=0;j<info->order;j+=b->dim)
  339. if(vorbis_book_decodev_set(b,lsp+j,&vd->opb,b->dim,-24)==-1)goto eop;
  340. for(j=0;j<info->order;){
  341. for(k=0;k<b->dim;k++,j++)lsp[j]+=last;
  342. last=lsp[j-1];
  343. }
  344. lsp[info->order]=amp;
  345. return(lsp);
  346. }
  347. }
  348. eop:
  349. return(NULL);
  350. }
  351. int floor0_inverse2(vorbis_dsp_state *vd,vorbis_info_floor *i,
  352. ogg_int32_t *lsp,ogg_int32_t *out){
  353. vorbis_info_floor0 *info=(vorbis_info_floor0 *)i;
  354. codec_setup_info *ci=(codec_setup_info *)vd->vi->codec_setup;
  355. if(lsp){
  356. ogg_int32_t amp=lsp[info->order];
  357. /* take the coefficients back to a spectral envelope curve */
  358. vorbis_lsp_to_curve(out,ci->blocksizes[vd->W]/2,info->barkmap,
  359. lsp,info->order,amp,info->ampdB,
  360. info->rate>>1);
  361. return(1);
  362. }
  363. memset(out,0,sizeof(*out)*ci->blocksizes[vd->W]/2);
  364. return(0);
  365. }