1 /*
2  ** Copyright 2003-2010, VisualOn, Inc.
3  **
4  ** Licensed under the Apache License, Version 2.0 (the "License");
5  ** you may not use this file except in compliance with the License.
6  ** You may obtain a copy of the License at
7  **
8  **     http://www.apache.org/licenses/LICENSE-2.0
9  **
10  ** Unless required by applicable law or agreed to in writing, software
11  ** distributed under the License is distributed on an "AS IS" BASIS,
12  ** WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13  ** See the License for the specific language governing permissions and
14  ** limitations under the License.
15  */
16 
17 /*****************************************************************************
18  *  This file contains operations in double precision.                       *
19  *  These operations are not standard double precision operations.           *
20  *  They are used where single precision is not enough but the full 32 bits  *
21  *  precision is not necessary. For example, the function Div_32() has a     *
22  *  24 bits precision which is enough for our purposes.                      *
23  *                                                                           *
24  *  The double precision numbers use a special representation:               *
25  *                                                                           *
26  *     L_32 = hi<<16 + lo<<1                                                 *
27  *                                                                           *
28  *  L_32 is a 32 bit integer.                                                *
29  *  hi and lo are 16 bit signed integers.                                    *
30  *  As the low part also contains the sign, this allows fast multiplication. *
31  *                                                                           *
32  *      0x8000 0000 <= L_32 <= 0x7fff fffe.                                  *
33  *                                                                           *
34  *  We will use DPF (Double Precision Format )in this file to specify        *
35  *  this special format.                                                     *
36  *****************************************************************************
37 */
38 #include "typedef.h"
39 #include "basic_op.h"
40 #include "oper_32b.h"
41 
42 /*****************************************************************************
43  *                                                                           *
44  *  Function L_Extract()                                                     *
45  *                                                                           *
46  *  Extract from a 32 bit integer two 16 bit DPF.                            *
47  *                                                                           *
48  *  Arguments:                                                               *
49  *                                                                           *
50  *   L_32      : 32 bit integer.                                             *
51  *               0x8000 0000 <= L_32 <= 0x7fff ffff.                         *
52  *   hi        : b16 to b31 of L_32                                          *
53  *   lo        : (L_32 - hi<<16)>>1                                          *
54  *****************************************************************************
55 */
56 
VO_L_Extract(Word32 L_32,Word16 * hi,Word16 * lo)57 __inline void VO_L_Extract (Word32 L_32, Word16 *hi, Word16 *lo)
58 {
59     *hi = (Word16)(L_32 >> 16);
60     *lo = (Word16)((L_32 & 0xffff) >> 1);
61     return;
62 }
63 
64 /*****************************************************************************
65  *                                                                           *
66  *  Function L_Comp()                                                        *
67  *                                                                           *
68  *  Compose from two 16 bit DPF a 32 bit integer.                            *
69  *                                                                           *
70  *     L_32 = hi<<16 + lo<<1                                                 *
71  *                                                                           *
72  *  Arguments:                                                               *
73  *                                                                           *
74  *   hi        msb                                                           *
75  *   lo        lsf (with sign)                                               *
76  *                                                                           *
77  *   Return Value :                                                          *
78  *                                                                           *
79  *             32 bit long signed integer (Word32) whose value falls in the  *
80  *             range : 0x8000 0000 <= L_32 <= 0x7fff fff0.                   *
81  *                                                                           *
82  *****************************************************************************
83 */
84 
L_Comp(Word16 hi,Word16 lo)85 Word32 L_Comp (Word16 hi, Word16 lo)
86 {
87     Word32 L_32;
88 
89     L_32 = L_deposit_h (hi);
90 
91     return (L_mac (L_32, lo, 1));       /* = hi<<16 + lo<<1 */
92 }
93 
94 /*****************************************************************************
95  * Function Mpy_32()                                                         *
96  *                                                                           *
97  *   Multiply two 32 bit integers (DPF). The result is divided by 2**31      *
98  *                                                                           *
99  *   L_32 = (hi1*hi2)<<1 + ( (hi1*lo2)>>15 + (lo1*hi2)>>15 )<<1              *
100  *                                                                           *
101  *   This operation can also be viewed as the multiplication of two Q31      *
102  *   number and the result is also in Q31.                                   *
103  *                                                                           *
104  * Arguments:                                                                *
105  *                                                                           *
106  *  hi1         hi part of first number                                      *
107  *  lo1         lo part of first number                                      *
108  *  hi2         hi part of second number                                     *
109  *  lo2         lo part of second number                                     *
110  *                                                                           *
111  *****************************************************************************
112 */
113 
Mpy_32(Word16 hi1,Word16 lo1,Word16 hi2,Word16 lo2)114 __inline Word32  Mpy_32 (Word16 hi1, Word16 lo1, Word16 hi2, Word16 lo2)
115 {
116     Word32 L_32;
117     L_32 = (hi1 * hi2);
118     L_32 += (hi1 * lo2) >> 15;
119     L_32 += (lo1 * hi2) >> 15;
120     L_32 <<= 1;
121 
122     return (L_32);
123 }
124 
125 /*****************************************************************************
126  * Function Mpy_32_16()                                                      *
127  *                                                                           *
128  *   Multiply a 16 bit integer by a 32 bit (DPF). The result is divided      *
129  *   by 2**15                                                                *
130  *                                                                           *
131  *                                                                           *
132  *   L_32 = (hi1*lo2)<<1 + ((lo1*lo2)>>15)<<1                                *
133  *                                                                           *
134  * Arguments:                                                                *
135  *                                                                           *
136  *  hi          hi part of 32 bit number.                                    *
137  *  lo          lo part of 32 bit number.                                    *
138  *  n           16 bit number.                                               *
139  *                                                                           *
140  *****************************************************************************
141 */
142 
Mpy_32_16(Word16 hi,Word16 lo,Word16 n)143 __inline Word32 Mpy_32_16 (Word16 hi, Word16 lo, Word16 n)
144 {
145     Word32 L_32;
146 
147     L_32 = (hi * n)<<1;
148     L_32 += (((lo * n)>>15)<<1);
149 
150     return (L_32);
151 }
152 
153 /*****************************************************************************
154  *                                                                           *
155  *   Function Name : Div_32                                                  *
156  *                                                                           *
157  *   Purpose :                                                               *
158  *             Fractional integer division of two 32 bit numbers.            *
159  *             L_num / L_denom.                                              *
160  *             L_num and L_denom must be positive and L_num < L_denom.       *
161  *             L_denom = denom_hi<<16 + denom_lo<<1                          *
162  *             denom_hi is a normalize number.                               *
163  *                                                                           *
164  *   Inputs :                                                                *
165  *                                                                           *
166  *    L_num                                                                  *
167  *             32 bit long signed integer (Word32) whose value falls in the  *
168  *             range : 0x0000 0000 < L_num < L_denom                         *
169  *                                                                           *
170  *    L_denom = denom_hi<<16 + denom_lo<<1      (DPF)                        *
171  *                                                                           *
172  *       denom_hi                                                            *
173  *             16 bit positive normalized integer whose value falls in the   *
174  *             range : 0x4000 < hi < 0x7fff                                  *
175  *       denom_lo                                                            *
176  *             16 bit positive integer whose value falls in the              *
177  *             range : 0 < lo < 0x7fff                                       *
178  *                                                                           *
179  *   Return Value :                                                          *
180  *                                                                           *
181  *    L_div                                                                  *
182  *             32 bit long signed integer (Word32) whose value falls in the  *
183  *             range : 0x0000 0000 <= L_div <= 0x7fff ffff.                  *
184  *                                                                           *
185  *  Algorithm:                                                               *
186  *                                                                           *
187  *  - find = 1/L_denom.                                                      *
188  *      First approximation: approx = 1 / denom_hi                           *
189  *      1/L_denom = approx * (2.0 - L_denom * approx )                       *
190  *                                                                           *
191  *  -  result = L_num * (1/L_denom)                                          *
192  *****************************************************************************
193 */
194 
Div_32(Word32 L_num,Word16 denom_hi,Word16 denom_lo)195 Word32 Div_32 (Word32 L_num, Word16 denom_hi, Word16 denom_lo)
196 {
197     Word16 approx, hi, lo, n_hi, n_lo;
198     Word32 L_32;
199 
200     /* First approximation: 1 / L_denom = 1/denom_hi */
201 
202     approx = div_s ((Word16) 0x3fff, denom_hi);
203 
204     /* 1/L_denom = approx * (2.0 - L_denom * approx) */
205 
206     L_32 = Mpy_32_16 (denom_hi, denom_lo, approx);
207 
208     L_32 = L_sub ((Word32) 0x7fffffffL, L_32);
209     hi = L_32 >> 16;
210     lo = (L_32 & 0xffff) >> 1;
211 
212     L_32 = Mpy_32_16 (hi, lo, approx);
213 
214     /* L_num * (1/L_denom) */
215     hi = L_32 >> 16;
216     lo = (L_32 & 0xffff) >> 1;
217     VO_L_Extract (L_num, &n_hi, &n_lo);
218     L_32 = Mpy_32 (n_hi, n_lo, hi, lo);
219     L_32 = L_shl2(L_32, 2);
220 
221     return (L_32);
222 }
223 
224