~ubuntu-branches/ubuntu/precise/linux-ti-omap4/precise

« back to all changes in this revision

Viewing changes to drivers/staging/brcm80211/util/qmath.c

  • Committer: Bazaar Package Importer
  • Author(s): Paolo Pisati
  • Date: 2011-06-29 15:23:51 UTC
  • mfrom: (26.1.1 natty-proposed)
  • Revision ID: james.westby@ubuntu.com-20110629152351-xs96tm303d95rpbk
Tags: 3.0.0-1200.2
* Rebased against 3.0.0-6.7
* BSP from TI based on 3.0.0

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
/*
2
 
 * Copyright (c) 2010 Broadcom Corporation
3
 
 *
4
 
 * Permission to use, copy, modify, and/or distribute this software for any
5
 
 * purpose with or without fee is hereby granted, provided that the above
6
 
 * copyright notice and this permission notice appear in all copies.
7
 
 *
8
 
 * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
9
 
 * WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
10
 
 * MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY
11
 
 * SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
12
 
 * WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN ACTION
13
 
 * OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF OR IN
14
 
 * CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
15
 
 */
16
 
 
17
 
#include <linux/types.h>
18
 
#include "qmath.h"
19
 
 
20
 
/*
21
 
Description: This function saturate input 32 bit number into a 16 bit number.
22
 
If input number is greater than 0x7fff then output is saturated to 0x7fff.
23
 
else if input number is less than 0xffff8000 then output is saturated to 0xffff8000
24
 
else output is same as input.
25
 
*/
26
 
s16 qm_sat32(s32 op)
27
 
{
28
 
        s16 result;
29
 
        if (op > (s32) 0x7fff) {
30
 
                result = 0x7fff;
31
 
        } else if (op < (s32) 0xffff8000) {
32
 
                result = (s16) (0x8000);
33
 
        } else {
34
 
                result = (s16) op;
35
 
        }
36
 
        return result;
37
 
}
38
 
 
39
 
/*
40
 
Description: This function multiply two input 16 bit numbers and return the 32 bit result.
41
 
This multiplication is similar to compiler multiplication. This operation is defined if
42
 
16 bit multiplication on the processor platform is cheaper than 32 bit multiplication (as
43
 
the most of qmath functions can be replaced with processor intrinsic instructions).
44
 
*/
45
 
s32 qm_mul321616(s16 op1, s16 op2)
46
 
{
47
 
        return (s32) (op1) * (s32) (op2);
48
 
}
49
 
 
50
 
/*
51
 
Description: This function make 16 bit multiplication and return the result in 16 bits.
52
 
To fit the result into 16 bits the 32 bit multiplication result is right
53
 
shifted by 16 bits.
54
 
*/
55
 
s16 qm_mul16(s16 op1, s16 op2)
56
 
{
57
 
        s32 result;
58
 
        result = ((s32) (op1) * (s32) (op2));
59
 
        return (s16) (result >> 16);
60
 
}
61
 
 
62
 
/*
63
 
Description: This function multiply two 16 bit numbers and return the result in 32 bits.
64
 
This function remove the extra sign bit created by the multiplication by leftshifting the
65
 
32 bit multiplication result by 1 bit before returning the result. So the output is
66
 
twice that of compiler multiplication. (i.e. qm_muls321616(2,3)=12).
67
 
When both input 16 bit numbers are 0x8000, then the result is saturated to 0x7fffffff.
68
 
*/
69
 
s32 qm_muls321616(s16 op1, s16 op2)
70
 
{
71
 
        s32 result;
72
 
        if (op1 == (s16) (0x8000) && op2 == (s16) (0x8000)) {
73
 
                result = 0x7fffffff;
74
 
        } else {
75
 
                result = ((s32) (op1) * (s32) (op2));
76
 
                result = result << 1;
77
 
        }
78
 
        return result;
79
 
}
80
 
 
81
 
/*
82
 
Description: This function make 16 bit unsigned multiplication. To fit the output into
83
 
16 bits the 32 bit multiplication result is right shifted by 16 bits.
84
 
*/
85
 
u16 qm_mulu16(u16 op1, u16 op2)
86
 
{
87
 
        return (u16) (((u32) op1 * (u32) op2) >> 16);
88
 
}
89
 
 
90
 
/*
91
 
Description: This function make 16 bit multiplication and return the result in 16 bits.
92
 
To fit the multiplication result into 16 bits the multiplication result is right shifted by
93
 
15 bits. Right shifting 15 bits instead of 16 bits is done to remove the extra sign bit formed
94
 
due to the multiplication.
95
 
When both the 16bit inputs are 0x8000 then the output is saturated to 0x7fffffff.
96
 
*/
97
 
s16 qm_muls16(s16 op1, s16 op2)
98
 
{
99
 
        s32 result;
100
 
        if (op1 == (s16) 0x8000 && op2 == (s16) 0x8000) {
101
 
                result = 0x7fffffff;
102
 
        } else {
103
 
                result = ((s32) (op1) * (s32) (op2));
104
 
        }
105
 
        return (s16) (result >> 15);
106
 
}
107
 
 
108
 
/*
109
 
Description: This function add two 32 bit numbers and return the 32bit result.
110
 
If the result overflow 32 bits, the output will be saturated to 32bits.
111
 
*/
112
 
s32 qm_add32(s32 op1, s32 op2)
113
 
{
114
 
        s32 result;
115
 
        result = op1 + op2;
116
 
        if (op1 < 0 && op2 < 0 && result > 0) {
117
 
                result = 0x80000000;
118
 
        } else if (op1 > 0 && op2 > 0 && result < 0) {
119
 
                result = 0x7fffffff;
120
 
        }
121
 
        return result;
122
 
}
123
 
 
124
 
/*
125
 
Description: This function add two 16 bit numbers and return the 16bit result.
126
 
If the result overflow 16 bits, the output will be saturated to 16bits.
127
 
*/
128
 
s16 qm_add16(s16 op1, s16 op2)
129
 
{
130
 
        s16 result;
131
 
        s32 temp = (s32) op1 + (s32) op2;
132
 
        if (temp > (s32) 0x7fff) {
133
 
                result = (s16) 0x7fff;
134
 
        } else if (temp < (s32) 0xffff8000) {
135
 
                result = (s16) 0xffff8000;
136
 
        } else {
137
 
                result = (s16) temp;
138
 
        }
139
 
        return result;
140
 
}
141
 
 
142
 
/*
143
 
Description: This function make 16 bit subtraction and return the 16bit result.
144
 
If the result overflow 16 bits, the output will be saturated to 16bits.
145
 
*/
146
 
s16 qm_sub16(s16 op1, s16 op2)
147
 
{
148
 
        s16 result;
149
 
        s32 temp = (s32) op1 - (s32) op2;
150
 
        if (temp > (s32) 0x7fff) {
151
 
                result = (s16) 0x7fff;
152
 
        } else if (temp < (s32) 0xffff8000) {
153
 
                result = (s16) 0xffff8000;
154
 
        } else {
155
 
                result = (s16) temp;
156
 
        }
157
 
        return result;
158
 
}
159
 
 
160
 
/*
161
 
Description: This function make 32 bit subtraction and return the 32bit result.
162
 
If the result overflow 32 bits, the output will be saturated to 32bits.
163
 
*/
164
 
s32 qm_sub32(s32 op1, s32 op2)
165
 
{
166
 
        s32 result;
167
 
        result = op1 - op2;
168
 
        if (op1 >= 0 && op2 < 0 && result < 0) {
169
 
                result = 0x7fffffff;
170
 
        } else if (op1 < 0 && op2 > 0 && result > 0) {
171
 
                result = 0x80000000;
172
 
        }
173
 
        return result;
174
 
}
175
 
 
176
 
/*
177
 
Description: This function multiply input 16 bit numbers and accumulate the result
178
 
into the input 32 bit number and return the 32 bit accumulated result.
179
 
If the accumulation result in overflow, then the output will be saturated.
180
 
*/
181
 
s32 qm_mac321616(s32 acc, s16 op1, s16 op2)
182
 
{
183
 
        s32 result;
184
 
        result = qm_add32(acc, qm_mul321616(op1, op2));
185
 
        return result;
186
 
}
187
 
 
188
 
/*
189
 
Description: This function make a 32 bit saturated left shift when the specified shift
190
 
is +ve. This function will make a 32 bit right shift when the specified shift is -ve.
191
 
This function return the result after shifting operation.
192
 
*/
193
 
s32 qm_shl32(s32 op, int shift)
194
 
{
195
 
        int i;
196
 
        s32 result;
197
 
        result = op;
198
 
        if (shift > 31)
199
 
                shift = 31;
200
 
        else if (shift < -31)
201
 
                shift = -31;
202
 
        if (shift >= 0) {
203
 
                for (i = 0; i < shift; i++) {
204
 
                        result = qm_add32(result, result);
205
 
                }
206
 
        } else {
207
 
                result = result >> (-shift);
208
 
        }
209
 
        return result;
210
 
}
211
 
 
212
 
/*
213
 
Description: This function make a 32 bit right shift when shift is +ve.
214
 
This function make a 32 bit saturated left shift when shift is -ve. This function
215
 
return the result of the shift operation.
216
 
*/
217
 
s32 qm_shr32(s32 op, int shift)
218
 
{
219
 
        return qm_shl32(op, -shift);
220
 
}
221
 
 
222
 
/*
223
 
Description: This function make a 16 bit saturated left shift when the specified shift
224
 
is +ve. This function will make a 16 bit right shift when the specified shift is -ve.
225
 
This function return the result after shifting operation.
226
 
*/
227
 
s16 qm_shl16(s16 op, int shift)
228
 
{
229
 
        int i;
230
 
        s16 result;
231
 
        result = op;
232
 
        if (shift > 15)
233
 
                shift = 15;
234
 
        else if (shift < -15)
235
 
                shift = -15;
236
 
        if (shift > 0) {
237
 
                for (i = 0; i < shift; i++) {
238
 
                        result = qm_add16(result, result);
239
 
                }
240
 
        } else {
241
 
                result = result >> (-shift);
242
 
        }
243
 
        return result;
244
 
}
245
 
 
246
 
/*
247
 
Description: This function make a 16 bit right shift when shift is +ve.
248
 
This function make a 16 bit saturated left shift when shift is -ve. This function
249
 
return the result of the shift operation.
250
 
*/
251
 
s16 qm_shr16(s16 op, int shift)
252
 
{
253
 
        return qm_shl16(op, -shift);
254
 
}
255
 
 
256
 
/*
257
 
Description: This function return the number of redundant sign bits in a 16 bit number.
258
 
Example: qm_norm16(0x0080) = 7.
259
 
*/
260
 
s16 qm_norm16(s16 op)
261
 
{
262
 
        u16 u16extraSignBits;
263
 
        if (op == 0) {
264
 
                return 15;
265
 
        } else {
266
 
                u16extraSignBits = 0;
267
 
                while ((op >> 15) == (op >> 14)) {
268
 
                        u16extraSignBits++;
269
 
                        op = op << 1;
270
 
                }
271
 
        }
272
 
        return u16extraSignBits;
273
 
}
274
 
 
275
 
/*
276
 
Description: This function return the number of redundant sign bits in a 32 bit number.
277
 
Example: qm_norm32(0x00000080) = 23
278
 
*/
279
 
s16 qm_norm32(s32 op)
280
 
{
281
 
        u16 u16extraSignBits;
282
 
        if (op == 0) {
283
 
                return 31;
284
 
        } else {
285
 
                u16extraSignBits = 0;
286
 
                while ((op >> 31) == (op >> 30)) {
287
 
                        u16extraSignBits++;
288
 
                        op = op << 1;
289
 
                }
290
 
        }
291
 
        return u16extraSignBits;
292
 
}
293
 
 
294
 
/*
295
 
Description: This function divide two 16 bit unsigned numbers.
296
 
The numerator should be less than denominator. So the quotient is always less than 1.
297
 
This function return the quotient in q.15 format.
298
 
*/
299
 
s16 qm_div_s(s16 num, s16 denom)
300
 
{
301
 
        s16 var_out;
302
 
        s16 iteration;
303
 
        s32 L_num;
304
 
        s32 L_denom;
305
 
        L_num = (num) << 15;
306
 
        L_denom = (denom) << 15;
307
 
        for (iteration = 0; iteration < 15; iteration++) {
308
 
                L_num <<= 1;
309
 
                if (L_num >= L_denom) {
310
 
                        L_num = qm_sub32(L_num, L_denom);
311
 
                        L_num = qm_add32(L_num, 1);
312
 
                }
313
 
        }
314
 
        var_out = (s16) (L_num & 0x7fff);
315
 
        return var_out;
316
 
}
317
 
 
318
 
/*
319
 
Description: This function compute the absolute value of a 16 bit number.
320
 
*/
321
 
s16 qm_abs16(s16 op)
322
 
{
323
 
        if (op < 0) {
324
 
                if (op == (s16) 0xffff8000) {
325
 
                        return 0x7fff;
326
 
                } else {
327
 
                        return -op;
328
 
                }
329
 
        } else {
330
 
                return op;
331
 
        }
332
 
}
333
 
 
334
 
/*
335
 
Description: This function divide two 16 bit numbers.
336
 
The quotient is returned through return value.
337
 
The qformat of the quotient is returned through the pointer (qQuotient) passed
338
 
to this function. The qformat of quotient is adjusted appropriately such that
339
 
the quotient occupies all 16 bits.
340
 
*/
341
 
s16 qm_div16(s16 num, s16 denom, s16 *qQuotient)
342
 
{
343
 
        s16 sign;
344
 
        s16 nNum, nDenom;
345
 
        sign = num ^ denom;
346
 
        num = qm_abs16(num);
347
 
        denom = qm_abs16(denom);
348
 
        nNum = qm_norm16(num);
349
 
        nDenom = qm_norm16(denom);
350
 
        num = qm_shl16(num, nNum - 1);
351
 
        denom = qm_shl16(denom, nDenom);
352
 
        *qQuotient = nNum - 1 - nDenom + 15;
353
 
        if (sign >= 0) {
354
 
                return qm_div_s(num, denom);
355
 
        } else {
356
 
                return -qm_div_s(num, denom);
357
 
        }
358
 
}
359
 
 
360
 
/*
361
 
Description: This function compute absolute value of a 32 bit number.
362
 
*/
363
 
s32 qm_abs32(s32 op)
364
 
{
365
 
        if (op < 0) {
366
 
                if (op == (s32) 0x80000000) {
367
 
                        return 0x7fffffff;
368
 
                } else {
369
 
                        return -op;
370
 
                }
371
 
        } else {
372
 
                return op;
373
 
        }
374
 
}
375
 
 
376
 
/*
377
 
Description: This function divide two 32 bit numbers. The division is performed
378
 
by considering only important 16 bits in 32 bit numbers.
379
 
The quotient is returned through return value.
380
 
The qformat of the quotient is returned through the pointer (qquotient) passed
381
 
to this function. The qformat of quotient is adjusted appropriately such that
382
 
the quotient occupies all 16 bits.
383
 
*/
384
 
s16 qm_div163232(s32 num, s32 denom, s16 *qquotient)
385
 
{
386
 
        s32 sign;
387
 
        s16 nNum, nDenom;
388
 
        sign = num ^ denom;
389
 
        num = qm_abs32(num);
390
 
        denom = qm_abs32(denom);
391
 
        nNum = qm_norm32(num);
392
 
        nDenom = qm_norm32(denom);
393
 
        num = qm_shl32(num, nNum - 1);
394
 
        denom = qm_shl32(denom, nDenom);
395
 
        *qquotient = nNum - 1 - nDenom + 15;
396
 
        if (sign >= 0) {
397
 
                return qm_div_s((s16) (num >> 16), (s16) (denom >> 16));
398
 
        } else {
399
 
                return -qm_div_s((s16) (num >> 16), (s16) (denom >> 16));
400
 
        }
401
 
}
402
 
 
403
 
/*
404
 
Description: This function multiply a 32 bit number with a 16 bit number.
405
 
The multiplicaton result is right shifted by 16 bits to fit the result
406
 
into 32 bit output.
407
 
*/
408
 
s32 qm_mul323216(s32 op1, s16 op2)
409
 
{
410
 
        s16 hi;
411
 
        u16 lo;
412
 
        s32 result;
413
 
        hi = op1 >> 16;
414
 
        lo = (s16) (op1 & 0xffff);
415
 
        result = qm_mul321616(hi, op2);
416
 
        result = result + (qm_mulsu321616(op2, lo) >> 16);
417
 
        return result;
418
 
}
419
 
 
420
 
/*
421
 
Description: This function multiply signed 16 bit number with unsigned 16 bit number and return
422
 
the result in 32 bits.
423
 
*/
424
 
s32 qm_mulsu321616(s16 op1, u16 op2)
425
 
{
426
 
        return (s32) (op1) * op2;
427
 
}
428
 
 
429
 
/*
430
 
Description: This function multiply 32 bit number with 16 bit number. The multiplication result is
431
 
right shifted by 15 bits to fit the result into 32 bits. Right shifting by only 15 bits instead of
432
 
16 bits is done to remove the extra sign bit formed by multiplication from the return value.
433
 
When the input numbers are 0x80000000, 0x8000 the return value is saturated to 0x7fffffff.
434
 
*/
435
 
s32 qm_muls323216(s32 op1, s16 op2)
436
 
{
437
 
        s16 hi;
438
 
        u16 lo;
439
 
        s32 result;
440
 
        hi = op1 >> 16;
441
 
        lo = (s16) (op1 & 0xffff);
442
 
        result = qm_muls321616(hi, op2);
443
 
        result = qm_add32(result, (qm_mulsu321616(op2, lo) >> 15));
444
 
        return result;
445
 
}
446
 
 
447
 
/*
448
 
Description: This function multiply two 32 bit numbers. The multiplication result is right
449
 
shifted by 32 bits to fit the multiplication result into 32 bits. The right shifted
450
 
multiplication result is returned as output.
451
 
*/
452
 
s32 qm_mul32(s32 a, s32 b)
453
 
{
454
 
        s16 hi1, hi2;
455
 
        u16 lo1, lo2;
456
 
        s32 result;
457
 
        hi1 = a >> 16;
458
 
        hi2 = b >> 16;
459
 
        lo1 = (u16) (a & 0xffff);
460
 
        lo2 = (u16) (b & 0xffff);
461
 
        result = qm_mul321616(hi1, hi2);
462
 
        result = result + (qm_mulsu321616(hi1, lo2) >> 16);
463
 
        result = result + (qm_mulsu321616(hi2, lo1) >> 16);
464
 
        return result;
465
 
}
466
 
 
467
 
/*
468
 
Description: This function multiply two 32 bit numbers. The multiplication result is
469
 
right shifted by 31 bits to fit the multiplication result into 32 bits. The right
470
 
shifted multiplication result is returned as output. Right shifting by only 31 bits
471
 
instead of 32 bits is done to remove the extra sign bit formed by multiplication.
472
 
When the input numbers are 0x80000000, 0x80000000 the return value is saturated to
473
 
0x7fffffff.
474
 
*/
475
 
s32 qm_muls32(s32 a, s32 b)
476
 
{
477
 
        s16 hi1, hi2;
478
 
        u16 lo1, lo2;
479
 
        s32 result;
480
 
        hi1 = a >> 16;
481
 
        hi2 = b >> 16;
482
 
        lo1 = (u16) (a & 0xffff);
483
 
        lo2 = (u16) (b & 0xffff);
484
 
        result = qm_muls321616(hi1, hi2);
485
 
        result = qm_add32(result, (qm_mulsu321616(hi1, lo2) >> 15));
486
 
        result = qm_add32(result, (qm_mulsu321616(hi2, lo1) >> 15));
487
 
        result = qm_add32(result, (qm_mulu16(lo1, lo2) >> 15));
488
 
        return result;
489
 
}
490
 
 
491
 
/* This table is log2(1+(i/32)) where i=[0:1:31], in q.15 format */
492
 
static const s16 log_table[] = {
493
 
        0,
494
 
        1455,
495
 
        2866,
496
 
        4236,
497
 
        5568,
498
 
        6863,
499
 
        8124,
500
 
        9352,
501
 
        10549,
502
 
        11716,
503
 
        12855,
504
 
        13968,
505
 
        15055,
506
 
        16117,
507
 
        17156,
508
 
        18173,
509
 
        19168,
510
 
        20143,
511
 
        21098,
512
 
        22034,
513
 
        22952,
514
 
        23852,
515
 
        24736,
516
 
        25604,
517
 
        26455,
518
 
        27292,
519
 
        28114,
520
 
        28922,
521
 
        29717,
522
 
        30498,
523
 
        31267,
524
 
        32024
525
 
};
526
 
 
527
 
#define LOG_TABLE_SIZE 32       /* log_table size */
528
 
#define LOG2_LOG_TABLE_SIZE 5   /* log2(log_table size) */
529
 
#define Q_LOG_TABLE 15          /* qformat of log_table */
530
 
#define LOG10_2         19728   /* log10(2) in q.16 */
531
 
 
532
 
/*
533
 
Description:
534
 
This routine takes the input number N and its q format qN and compute
535
 
the log10(N). This routine first normalizes the input no N.     Then N is in mag*(2^x) format.
536
 
mag is any number in the range 2^30-(2^31 - 1). Then log2(mag * 2^x) = log2(mag) + x is computed.
537
 
From that log10(mag * 2^x) = log2(mag * 2^x) * log10(2) is computed.
538
 
This routine looks the log2 value in the table considering LOG2_LOG_TABLE_SIZE+1 MSBs.
539
 
As the MSB is always 1, only next LOG2_OF_LOG_TABLE_SIZE MSBs are used for table lookup.
540
 
Next 16 MSBs are used for interpolation.
541
 
Inputs:
542
 
N - number to which log10 has to be found.
543
 
qN - q format of N
544
 
log10N - address where log10(N) will be written.
545
 
qLog10N - address where log10N qformat will be written.
546
 
Note/Problem:
547
 
For accurate results input should be in normalized or near normalized form.
548
 
*/
549
 
void qm_log10(s32 N, s16 qN, s16 *log10N, s16 *qLog10N)
550
 
{
551
 
        s16 s16norm, s16tableIndex, s16errorApproximation;
552
 
        u16 u16offset;
553
 
        s32 s32log;
554
 
 
555
 
        /* Logerithm of negative values is undefined.
556
 
         * assert N is greater than 0.
557
 
         */
558
 
        /* ASSERT(N > 0); */
559
 
 
560
 
        /* normalize the N. */
561
 
        s16norm = qm_norm32(N);
562
 
        N = N << s16norm;
563
 
 
564
 
        /* The qformat of N after normalization.
565
 
         * -30 is added to treat the no as between 1.0 to 2.0
566
 
         * i.e. after adding the -30 to the qformat the decimal point will be
567
 
         * just rigtht of the MSB. (i.e. after sign bit and 1st MSB). i.e.
568
 
         * at the right side of 30th bit.
569
 
         */
570
 
        qN = qN + s16norm - 30;
571
 
 
572
 
        /* take the table index as the LOG2_OF_LOG_TABLE_SIZE bits right of the MSB */
573
 
        s16tableIndex = (s16) (N >> (32 - (2 + LOG2_LOG_TABLE_SIZE)));
574
 
 
575
 
        /* remove the MSB. the MSB is always 1 after normalization. */
576
 
        s16tableIndex =
577
 
            s16tableIndex & (s16) ((1 << LOG2_LOG_TABLE_SIZE) - 1);
578
 
 
579
 
        /* remove the (1+LOG2_OF_LOG_TABLE_SIZE) MSBs in the N. */
580
 
        N = N & ((1 << (32 - (2 + LOG2_LOG_TABLE_SIZE))) - 1);
581
 
 
582
 
        /* take the offset as the 16 MSBS after table index.
583
 
         */
584
 
        u16offset = (u16) (N >> (32 - (2 + LOG2_LOG_TABLE_SIZE + 16)));
585
 
 
586
 
        /* look the log value in the table. */
587
 
        s32log = log_table[s16tableIndex];      /* q.15 format */
588
 
 
589
 
        /* interpolate using the offset. */
590
 
        s16errorApproximation = (s16) qm_mulu16(u16offset, (u16) (log_table[s16tableIndex + 1] - log_table[s16tableIndex]));    /* q.15 */
591
 
 
592
 
        s32log = qm_add16((s16) s32log, s16errorApproximation); /* q.15 format */
593
 
 
594
 
        /* adjust for the qformat of the N as
595
 
         * log2(mag * 2^x) = log2(mag) + x
596
 
         */
597
 
        s32log = qm_add32(s32log, ((s32) -qN) << 15);   /* q.15 format */
598
 
 
599
 
        /* normalize the result. */
600
 
        s16norm = qm_norm32(s32log);
601
 
 
602
 
        /* bring all the important bits into lower 16 bits */
603
 
        s32log = qm_shl32(s32log, s16norm - 16);        /* q.15+s16norm-16 format */
604
 
 
605
 
        /* compute the log10(N) by multiplying log2(N) with log10(2).
606
 
         * as log10(mag * 2^x) = log2(mag * 2^x) * log10(2)
607
 
         * log10N in q.15+s16norm-16+1 (LOG10_2 is in q.16)
608
 
         */
609
 
        *log10N = qm_muls16((s16) s32log, (s16) LOG10_2);
610
 
 
611
 
        /* write the q format of the result. */
612
 
        *qLog10N = 15 + s16norm - 16 + 1;
613
 
 
614
 
        return;
615
 
}
616
 
 
617
 
/*
618
 
Description:
619
 
This routine compute 1/N.
620
 
This routine reformates the given no N as N * 2^qN where N is in between 0.5 and 1.0
621
 
in q.15 format in 16 bits. So the problem now boils down to finding the inverse of a
622
 
q.15 no in 16 bits which is in the range of 0.5 to 1.0. The output is always between
623
 
2.0 to 1. So the output is 2.0 to 1.0 in q.30 format. Once the final output format is found
624
 
by taking the qN into account. Inverse is found with newton rapson method. Initially
625
 
inverse (x) is guessed as 1/0.75 (with appropriate sign). The new guess is calculated
626
 
using the formula x' = 2*x - N*x*x. After 4 or 5 iterations the inverse is very close to
627
 
inverse of N.
628
 
Inputs:
629
 
N - number to which 1/N has to be found.
630
 
qn - q format of N.
631
 
sqrtN - address where 1/N has to be written.
632
 
qsqrtN - address where q format of 1/N has to be written.
633
 
*/
634
 
#define qx 29
635
 
void qm_1byN(s32 N, s16 qN, s32 *result, s16 *qResult)
636
 
{
637
 
        s16 normN;
638
 
        s32 s32firstTerm, s32secondTerm, x;
639
 
        int i;
640
 
 
641
 
        normN = qm_norm32(N);
642
 
 
643
 
        /* limit N to least significant 16 bits. 15th bit is the sign bit. */
644
 
        N = qm_shl32(N, normN - 16);
645
 
        qN = qN + normN - 16 - 15;
646
 
        /* -15 is added to treat N as 16 bit q.15 number in the range from 0.5 to 1 */
647
 
 
648
 
        /* Take the initial guess as 1/0.75 in qx format with appropriate sign. */
649
 
        if (N >= 0) {
650
 
                x = (s32) ((1 / 0.75) * (1 << qx));
651
 
                /* input no is in the range 0.5 to 1. So 1/0.75 is taken as initial guess. */
652
 
        } else {
653
 
                x = (s32) ((1 / -0.75) * (1 << qx));
654
 
                /* input no is in the range -0.5 to -1. So 1/-0.75 is taken as initial guess. */
655
 
        }
656
 
 
657
 
        /* iterate the equation x = 2*x - N*x*x for 4 times. */
658
 
        for (i = 0; i < 4; i++) {
659
 
                s32firstTerm = qm_shl32(x, 1);  /* s32firstTerm = 2*x in q.29 */
660
 
                s32secondTerm =
661
 
                    qm_muls321616((s16) (s32firstTerm >> 16),
662
 
                                  (s16) (s32firstTerm >> 16));
663
 
                /* s32secondTerm = x*x in q.(29+1-16)*2+1 */
664
 
                s32secondTerm =
665
 
                    qm_muls321616((s16) (s32secondTerm >> 16), (s16) N);
666
 
                /* s32secondTerm = N*x*x in q.((29+1-16)*2+1)-16+15+1 i.e. in q.29 */
667
 
                x = qm_sub32(s32firstTerm, s32secondTerm);
668
 
                /* can be added directly as both are in q.29 */
669
 
        }
670
 
 
671
 
        /* Bring the x to q.30 format. */
672
 
        *result = qm_shl32(x, 1);
673
 
        /* giving the output in q.30 format for q.15 input in 16 bits. */
674
 
 
675
 
        /* compute the final q format of the result. */
676
 
        *qResult = -qN + 30;    /* adjusting the q format of actual output */
677
 
 
678
 
        return;
679
 
}
680
 
 
681
 
#undef qx