1e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard/*
2e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard ** Copyright 2003-2010, VisualOn, Inc.
3e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard **
4e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard ** Licensed under the Apache License, Version 2.0 (the "License");
5e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard ** you may not use this file except in compliance with the License.
6e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard ** You may obtain a copy of the License at
7e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard **
8e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard **     http://www.apache.org/licenses/LICENSE-2.0
9e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard **
10e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard ** Unless required by applicable law or agreed to in writing, software
11e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard ** distributed under the License is distributed on an "AS IS" BASIS,
12e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard ** WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard ** See the License for the specific language governing permissions and
14e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard ** limitations under the License.
15e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard */
16e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard/*******************************************************************************
17e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard	File:		oper_32b.c
18e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
19956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong	Content:	  This file contains operations in double precision.
20b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard
21956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong*******************************************************************************/
22956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
23956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong#include "typedef.h"
24956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong#include "basic_op.h"
25956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong#include "oper_32b.h"
26956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
27956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong/*****************************************************************************
28956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
29956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  Function L_Extract()                                                     *
30956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
31956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  Extract from a 32 bit integer two 16 bit DPF.                            *
32956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
33956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  Arguments:                                                               *
34956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
35956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   L_32      : 32 bit integer.                                             *
36956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *               0x8000 0000 <= L_32 <= 0x7fff ffff.                         *
37956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   hi        : b16 to b31 of L_32                                          *
38956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   lo        : (L_32 - hi<<16)>>1                                          *
39956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *****************************************************************************
40956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong*/
41956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
42956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dongvoid L_Extract (Word32 L_32, Word16 *hi, Word16 *lo)
43956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong{
44956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    *hi = extract_h (L_32);
45956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    *lo = extract_l (L_msu (L_shr (L_32, 1), *hi, 16384));
46956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    return;
47956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong}
48956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
49956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong/*****************************************************************************
50956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
51956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  Function L_Comp()                                                        *
52956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
53956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  Compose from two 16 bit DPF a 32 bit integer.                            *
54956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
55956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *     L_32 = hi<<16 + lo<<1                                                 *
56956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
57956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  Arguments:                                                               *
58956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
59956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   hi        msb                                                           *
60956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   lo        lsf (with sign)                                               *
61956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
62956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   Return Value :                                                          *
63956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
64956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             32 bit long signed integer (Word32) whose value falls in the  *
65956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             range : 0x8000 0000 <= L_32 <= 0x7fff fff0.                   *
66956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
67956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *****************************************************************************
68956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong*/
69956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
70956c553ab0ce72f8074ad0fda2ffd66a0305700cJames DongWord32 L_Comp (Word16 hi, Word16 lo)
71956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong{
72956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    Word32 L_32;
73956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
74956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    L_32 = L_deposit_h (hi);
75956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    return (L_mac (L_32, lo, 1));       /* = hi<<16 + lo<<1 */
76956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong}
77956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
78956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong/*****************************************************************************
79956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong * Function Mpy_32()                                                         *
80956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
81956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   Multiply two 32 bit integers (DPF). The result is divided by 2**31      *
82956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
83956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   L_32 = (hi1*hi2)<<1 + ( (hi1*lo2)>>15 + (lo1*hi2)>>15 )<<1              *
84956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
85956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   This operation can also be viewed as the multiplication of two Q31      *
86956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   number and the result is also in Q31.                                   *
87956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
88956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong * Arguments:                                                                *
89956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
90956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  hi1         hi part of first number                                      *
91956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  lo1         lo part of first number                                      *
92956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  hi2         hi part of second number                                     *
93956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  lo2         lo part of second number                                     *
94956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
95956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *****************************************************************************
96956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong*/
97956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
98956c553ab0ce72f8074ad0fda2ffd66a0305700cJames DongWord32 Mpy_32 (Word16 hi1, Word16 lo1, Word16 hi2, Word16 lo2)
99956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong{
100956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    Word32 L_32;
101956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
102956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    L_32 = L_mult (hi1, hi2);
103956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    L_32 = L_mac (L_32, mult (hi1, lo2), 1);
104956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    L_32 = L_mac (L_32, mult (lo1, hi2), 1);
105956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
106956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    return (L_32);
107956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong}
108956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
109956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong/*****************************************************************************
110956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong * Function Mpy_32_16()                                                      *
111956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
112956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   Multiply a 16 bit integer by a 32 bit (DPF). The result is divided      *
113956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   by 2**15                                                                *
114956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
115956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
116956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   L_32 = (hi1*lo2)<<1 + ((lo1*lo2)>>15)<<1                                *
117956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
118956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong * Arguments:                                                                *
119956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
120956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  hi          hi part of 32 bit number.                                    *
121956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  lo          lo part of 32 bit number.                                    *
122956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  n           16 bit number.                                               *
123956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
124956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *****************************************************************************
125956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong*/
126956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
127956c553ab0ce72f8074ad0fda2ffd66a0305700cJames DongWord32 Mpy_32_16 (Word16 hi, Word16 lo, Word16 n)
128956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong{
129956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    Word32 L_32;
130956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
131956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    L_32 = L_mult (hi, n);
132956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    L_32 = L_mac (L_32, mult (lo, n), 1);
133956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
134956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    return (L_32);
135956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong}
136956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
137956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong/*****************************************************************************
138956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
139956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   Function Name : Div_32                                                  *
140956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
141956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   Purpose :                                                               *
142956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             Fractional integer division of two 32 bit numbers.            *
143956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             L_num / L_denom.                                              *
144956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             L_num and L_denom must be positive and L_num < L_denom.       *
145956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             L_denom = denom_hi<<16 + denom_lo<<1                          *
146956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             denom_hi is a normalize number.                               *
147956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
148956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   Inputs :                                                                *
149956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
150956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *    L_num                                                                  *
151956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             32 bit long signed integer (Word32) whose value falls in the  *
152956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             range : 0x0000 0000 < L_num < L_denom                         *
153956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
154956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *    L_denom = denom_hi<<16 + denom_lo<<1      (DPF)                        *
155956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
156956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *       denom_hi                                                            *
157956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             16 bit positive normalized integer whose value falls in the   *
158956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             range : 0x4000 < hi < 0x7fff                                  *
159956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *       denom_lo                                                            *
160956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             16 bit positive integer whose value falls in the              *
161956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             range : 0 < lo < 0x7fff                                       *
162956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
163956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *   Return Value :                                                          *
164956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
165956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *    L_div                                                                  *
166956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             32 bit long signed integer (Word32) whose value falls in the  *
167956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *             range : 0x0000 0000 <= L_div <= 0x7fff ffff.                  *
168956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
169956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  Algorithm:                                                               *
170956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
171956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  - find = 1/L_denom.                                                      *
172956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *      First approximation: approx = 1 / denom_hi                           *
173956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *      1/L_denom = approx * (2.0 - L_denom * approx )                       *
174956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *                                                                           *
175956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *  -  result = L_num * (1/L_denom)                                          *
176956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong *****************************************************************************
177956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong*/
178956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
179956c553ab0ce72f8074ad0fda2ffd66a0305700cJames DongWord32 Div_32 (Word32 L_num, Word32 denom)
180956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong{
181956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    Word16 approx;
182956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    Word32 L_32;
183956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    /* First approximation: 1 / L_denom = 1/denom_hi */
184956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
185956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    approx = div_s ((Word16) 0x3fff, denom >> 16);
186956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
187956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    /* 1/L_denom = approx * (2.0 - L_denom * approx) */
188956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
189956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    L_32 = L_mpy_ls (denom, approx);
190956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
191956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    L_32 = L_sub ((Word32) 0x7fffffffL, L_32);
192956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
193956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong	L_32 = L_mpy_ls (L_32, approx);
194956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    /* L_num * (1/L_denom) */
195956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
196956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong	L_32 = MULHIGH(L_32, L_num);
197956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    L_32 = L_shl (L_32, 3);
198956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong
199956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong    return (L_32);
200956c553ab0ce72f8074ad0fda2ffd66a0305700cJames Dong}
201e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
202e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard/*!
203b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard
204b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard  \brief  calculates the log dualis times 4 of argument
205e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard          iLog4(x) = (Word32)(4 * log(value)/log(2.0))
206e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
207e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  \return ilog4 value
208b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard
209e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard*/
210e2e838afcf03e603a41a0455846eaf9614537c16Mans RullgardWord16 iLog4(Word32 value)
211e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard{
212e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  Word16 iLog4;
213e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
214e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  if(value != 0){
215e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    Word32 tmp;
216e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    Word16 tmp16;
217e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    iLog4 = norm_l(value);
218e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    tmp = (value << iLog4);
219e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    tmp16 = round16(tmp);
220e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    tmp = L_mult(tmp16, tmp16);
221e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    tmp16 = round16(tmp);
222e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    tmp = L_mult(tmp16, tmp16);
223e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    tmp16 = round16(tmp);
224e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
225e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    iLog4 = (-(iLog4 << 2) - norm_s(tmp16)) - 1;
226e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  }
227e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  else {
228b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard    iLog4 = -128; /* -(INT_BITS*4); */
229e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  }
230e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
231e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  return iLog4;
232e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard}
233e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
234e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard#define step(shift) \
235e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    if ((0x40000000l >> shift) + root <= value)       \
236e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    {                                                 \
237e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard        value -= (0x40000000l >> shift) + root;       \
238e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard        root = (root >> 1) | (0x40000000l >> shift);  \
239e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    } else {                                          \
240e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard        root = root >> 1;                             \
241e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    }
242e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
243e2e838afcf03e603a41a0455846eaf9614537c16Mans RullgardWord32 rsqrt(Word32 value,     /*!< Operand to square root (0.0 ... 1) */
244e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard             Word32 accuracy)  /*!< Number of valid bits that will be calculated */
245e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard{
246e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    Word32 root = 0;
247e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard	Word32 scale;
248e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
249e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard	if(value < 0)
250e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard		return 0;
251e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
252e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard	scale = norm_l(value);
253e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard	if(scale & 1) scale--;
254e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
255e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard	value <<= scale;
256e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
257e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard	step( 0); step( 2); step( 4); step( 6);
258e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    step( 8); step(10); step(12); step(14);
259e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    step(16); step(18); step(20); step(22);
260e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    step(24); step(26); step(28); step(30);
261e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
262e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    scale >>= 1;
263e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard	if (root < value)
264e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard        ++root;
265e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
266e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard	root >>= scale;
267e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard    return root* 46334;
268e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard}
269e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
270e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgardstatic const Word32 pow2Table[POW2_TABLE_SIZE] = {
271b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x7fffffff, 0x7fa765ad, 0x7f4f08ae, 0x7ef6e8da,
272b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x7e9f0606, 0x7e476009, 0x7deff6b6, 0x7d98c9e6,
273b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x7d41d96e, 0x7ceb2523, 0x7c94acde, 0x7c3e7073,
274b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x7be86fb9, 0x7b92aa88, 0x7b3d20b6, 0x7ae7d21a,
275b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x7a92be8b, 0x7a3de5df, 0x79e947ef, 0x7994e492,
276b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x7940bb9e, 0x78ecccec, 0x78991854, 0x78459dac,
277b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x77f25cce, 0x779f5591, 0x774c87cc, 0x76f9f359,
278b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x76a7980f, 0x765575c8, 0x76038c5b, 0x75b1dba2,
279b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x75606374, 0x750f23ab, 0x74be1c20, 0x746d4cac,
280b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x741cb528, 0x73cc556d, 0x737c2d55, 0x732c3cba,
281b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x72dc8374, 0x728d015d, 0x723db650, 0x71eea226,
282b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x719fc4b9, 0x71511de4, 0x7102ad80, 0x70b47368,
283b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x70666f76, 0x7018a185, 0x6fcb096f, 0x6f7da710,
284b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x6f307a41, 0x6ee382de, 0x6e96c0c3, 0x6e4a33c9,
285b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x6dfddbcc, 0x6db1b8a8, 0x6d65ca38, 0x6d1a1057,
286b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x6cce8ae1, 0x6c8339b2, 0x6c381ca6, 0x6bed3398,
287b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x6ba27e66, 0x6b57fce9, 0x6b0daeff, 0x6ac39485,
288b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x6a79ad56, 0x6a2ff94f, 0x69e6784d, 0x699d2a2c,
289b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x69540ec9, 0x690b2601, 0x68c26fb1, 0x6879ebb6,
290b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x683199ed, 0x67e97a34, 0x67a18c68, 0x6759d065,
291b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x6712460b, 0x66caed35, 0x6683c5c3, 0x663ccf92,
292b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x65f60a80, 0x65af766a, 0x6569132f, 0x6522e0ad,
293b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x64dcdec3, 0x64970d4f, 0x64516c2e, 0x640bfb41,
294b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x63c6ba64, 0x6381a978, 0x633cc85b, 0x62f816eb,
295b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x62b39509, 0x626f4292, 0x622b1f66, 0x61e72b65,
296b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x61a3666d, 0x615fd05f, 0x611c6919, 0x60d9307b,
297b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x60962665, 0x60534ab7, 0x60109d51, 0x5fce1e12,
298b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x5f8bccdb, 0x5f49a98c, 0x5f07b405, 0x5ec5ec26,
299b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x5e8451d0, 0x5e42e4e3, 0x5e01a540, 0x5dc092c7,
300b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x5d7fad59, 0x5d3ef4d7, 0x5cfe6923, 0x5cbe0a1c,
301b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x5c7dd7a4, 0x5c3dd19c, 0x5bfdf7e5, 0x5bbe4a61,
302b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x5b7ec8f2, 0x5b3f7377, 0x5b0049d4, 0x5ac14bea,
303b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x5a82799a, 0x5a43d2c6, 0x5a055751, 0x59c7071c,
304b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x5988e209, 0x594ae7fb, 0x590d18d3, 0x58cf7474,
305b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x5891fac1, 0x5854ab9b, 0x581786e6, 0x57da8c83,
306b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x579dbc57, 0x57611642, 0x57249a29, 0x56e847ef,
307b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x56ac1f75, 0x567020a0, 0x56344b52, 0x55f89f70,
308b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x55bd1cdb, 0x5581c378, 0x55469329, 0x550b8bd4,
309b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x54d0ad5b, 0x5495f7a1, 0x545b6a8b, 0x542105fd,
310b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x53e6c9db, 0x53acb607, 0x5372ca68, 0x533906e0,
311b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x52ff6b55, 0x52c5f7aa, 0x528cabc3, 0x52538786,
312b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x521a8ad7, 0x51e1b59a, 0x51a907b4, 0x5170810b,
313b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x51382182, 0x50ffe8fe, 0x50c7d765, 0x508fec9c,
314b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x50582888, 0x50208b0e, 0x4fe91413, 0x4fb1c37c,
315b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x4f7a9930, 0x4f439514, 0x4f0cb70c, 0x4ed5ff00,
316b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x4e9f6cd4, 0x4e69006e, 0x4e32b9b4, 0x4dfc988c,
317b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x4dc69cdd, 0x4d90c68b, 0x4d5b157e, 0x4d25899c,
318b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x4cf022ca, 0x4cbae0ef, 0x4c85c3f1, 0x4c50cbb8,
319b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x4c1bf829, 0x4be7492b, 0x4bb2bea5, 0x4b7e587d,
320b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x4b4a169c, 0x4b15f8e6, 0x4ae1ff43, 0x4aae299b,
321b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x4a7a77d5, 0x4a46e9d6, 0x4a137f88, 0x49e038d0,
322b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x49ad1598, 0x497a15c4, 0x4947393f, 0x49147fee,
323b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x48e1e9ba, 0x48af768a, 0x487d2646, 0x484af8d6,
324b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x4818ee22, 0x47e70611, 0x47b5408c, 0x47839d7b,
325b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x47521cc6, 0x4720be55, 0x46ef8210, 0x46be67e0,
326b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x468d6fae, 0x465c9961, 0x462be4e2, 0x45fb521a,
327b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x45cae0f2, 0x459a9152, 0x456a6323, 0x453a564d,
328b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x450a6abb, 0x44daa054, 0x44aaf702, 0x447b6ead,
329b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x444c0740, 0x441cc0a3, 0x43ed9ac0, 0x43be9580,
330b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x438fb0cb, 0x4360ec8d, 0x433248ae, 0x4303c517,
331b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x42d561b4, 0x42a71e6c, 0x4278fb2b, 0x424af7da,
332b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x421d1462, 0x41ef50ae, 0x41c1aca8, 0x41942839,
333b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x4166c34c, 0x41397dcc, 0x410c57a2, 0x40df50b8,
334b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard0x40b268fa, 0x4085a051, 0x4058f6a8, 0x402c6be9
335e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard};
336e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
337e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard/*!
338b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard
339b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard  \brief calculates 2 ^ (x/y) for x<=0, y > 0, x <= 32768 * y
340b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard
341e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  avoids integer division
342b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard
343b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard  \return
344e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard*/
345e2e838afcf03e603a41a0455846eaf9614537c16Mans RullgardWord32 pow2_xy(Word32 x, Word32 y)
346e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard{
347f8321d624f4bbbfcf01d59f346d3eb390f75f24cMartin Storsjo  UWord32 iPart;
348f8321d624f4bbbfcf01d59f346d3eb390f75f24cMartin Storsjo  UWord32 fPart;
349e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  Word32 res;
350e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  Word32 tmp, tmp2;
351e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  Word32 shift, shift2;
352e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
353e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  tmp2 = -x;
354e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  iPart = tmp2 / y;
355e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  fPart = tmp2 - iPart*y;
356e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  iPart = min(iPart,INT_BITS-1);
357e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard
358b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard  res = pow2Table[(POW2_TABLE_SIZE*fPart)/y] >> iPart;
359b676a05348e4c516fa8b57e33b10548e6142c3f8Mans Rullgard
360e2e838afcf03e603a41a0455846eaf9614537c16Mans Rullgard  return(res);
361891abc0ee089f2ba5b92dcc014e5efc2ef07f01eMartin Storsjo}
362