14f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber/* ------------------------------------------------------------------
24f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber * Copyright (C) 1998-2009 PacketVideo
34f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber *
44f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber * Licensed under the Apache License, Version 2.0 (the "License");
54f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber * you may not use this file except in compliance with the License.
64f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber * You may obtain a copy of the License at
74f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber *
84f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber *      http://www.apache.org/licenses/LICENSE-2.0
94f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber *
104f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber * Unless required by applicable law or agreed to in writing, software
114f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber * distributed under the License is distributed on an "AS IS" BASIS,
124f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either
134f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber * express or implied.
144f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber * See the License for the specific language governing permissions
154f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber * and limitations under the License.
164f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber * -------------------------------------------------------------------
174f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber */
184f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber/****************************************************************************************
194f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas HuberPortions of this file are derived from the following 3GPP standard:
204f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
214f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    3GPP TS 26.073
224f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    ANSI-C code for the Adaptive Multi-Rate (AMR) speech codec
234f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    Available from http://www.3gpp.org
244f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
254f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber(C) 2004, 3GPP Organizational Partners (ARIB, ATIS, CCSA, ETSI, TTA, TTC)
264f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas HuberPermission to distribute, modify and use this file under the standard license
274f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huberterms listed above has been obtained from the copyright holder.
284f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber****************************************************************************************/
294f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber/*
304f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
314f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber Pathname: ./audio/gsm-amr/c/src/sqrt_l.c
324f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
334f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber------------------------------------------------------------------------------
344f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber REVISION HISTORY
354f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
364f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
374f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber Description: Updated template. Changed function interface to pass in a
384f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber              pointer to overflow flag into the function instead of using a
394f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber              global flag. Changed name of an input pointer from "exp" to "pExp"
404f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber              for clarity. Removed inclusion of unwanted header files.
414f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
424f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber Description: Removed inclusion of sqrt_l.tab file. Changed the array name
434f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber              "table" to "sqrt_l_tbl". Fixed typos.
444f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
454f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber Who:                           Date:
464f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber Description:
474f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
484f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber------------------------------------------------------------------------------
494f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber*/
504f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
514f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber/*----------------------------------------------------------------------------
524f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber; INCLUDES
534f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber----------------------------------------------------------------------------*/
544f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber#include    "sqrt_l.h"
554f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber#include    "typedef.h"
564f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber#include    "basic_op.h"
574f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
584f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber/*----------------------------------------------------------------------------
594f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber; MACROS
604f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber; Define module specific macros here
614f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber----------------------------------------------------------------------------*/
624f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
634f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
644f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber/*----------------------------------------------------------------------------
654f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber; DEFINES
664f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber; Include all pre-processor statements here. Include conditional
674f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber; compile variables also.
684f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber----------------------------------------------------------------------------*/
694f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
704f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
714f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber/*----------------------------------------------------------------------------
724f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber; LOCAL FUNCTION DEFINITIONS
734f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber; Function Prototype declaration
744f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber----------------------------------------------------------------------------*/
754f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
764f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
774f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber/*----------------------------------------------------------------------------
784f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber; LOCAL STORE/BUFFER/POINTER DEFINITIONS
794f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber; Variable declaration - defined here and used outside this module
804f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber----------------------------------------------------------------------------*/
814f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
824f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
834f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber/*
844f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber------------------------------------------------------------------------------
854f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber FUNCTION NAME: sqrt_l_exp
864f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber------------------------------------------------------------------------------
874f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber INPUT AND OUTPUT DEFINITIONS
884f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
894f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber Inputs:
904f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    L_x = input value (Word32)
914f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    pExp = pointer to right shift to be applied to result
924f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    pOverflow = pointer to overflow flag
934f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
944f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber Outputs:
954f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    pOverflow -> if the Inv_sqrt operation resulted in an overflow.
964f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
974f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber Returns:
984f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    L_y = squareroot of L_x (Word32)
994f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1004f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber Global Variables Used:
1014f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    None.
1024f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1034f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber Local Variables Needed:
1044f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    None.
1054f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1064f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber------------------------------------------------------------------------------
1074f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber FUNCTION DESCRIPTION
1084f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1094f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber This function computes sqrt(L_x),  where  L_x is positive.
1104f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber If L_var is negative or zero, the result is 0
1114f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1124f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber The function sqrt(L_x) is approximated by a table and linear
1134f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber interpolation. The square root is computed using the
1144f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber following steps:
1154f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    1- Normalization of L_x.
1164f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    2- If exponent is even then shift right once.
1174f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    3- exponent = exponent/2
1184f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    4- i = bit25-b31 of L_x;  16<=i<=63  because of normalization.
1194f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    5- a = bit10-b24
1204f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    6- i -=16
1214f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    7- L_y = table[i]<<16 - (table[i] - table[i+1]) * a * 2
1224f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    8- return L_y and exponent so caller can do denormalization
1234f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1244f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber------------------------------------------------------------------------------
1254f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber REQUIREMENTS
1264f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1274f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber None.
1284f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1294f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber------------------------------------------------------------------------------
1304f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber REFERENCES
1314f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1324f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber sqrt_l.c, UMTS GSM AMR speech codec, R99 - Version 3.2.0, March 2, 2001
1334f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1344f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber------------------------------------------------------------------------------
1354f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber PSEUDO-CODE
1364f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1374f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas HuberWord32 sqrt_l_exp (     // o : output value
1384f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    Word32 L_x,         // i : input value
1394f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    Word16 *exp         // o : right shift to be applied to result
1404f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber)
1414f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber{
1424f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1434f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber//          y = sqrt(x)
1444f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber//          x = f * 2^-e,   0.5 <= f < 1   (normalization)
1454f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber//          y = sqrt(f) * 2^(-e/2)
1464f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber//
1474f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber//          a) e = 2k   --> y = sqrt(f)   * 2^-k  (k = e div 2,
1484f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber//                                                 0.707 <= sqrt(f) < 1)
1494f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber//          b) e = 2k+1 --> y = sqrt(f/2) * 2^-k  (k = e div 2,
1504f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber                                                 0.5 <= sqrt(f/2) < 0.707)
1514f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1524f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1534f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    Word16 e, i, a, tmp;
1544f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    Word32 L_y;
1554f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1564f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    if (L_x <= (Word32) 0)
1574f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    {
1584f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber        *exp = 0;
1594f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber        return (Word32) 0;
1604f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    }
1614f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1624f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber* The reference ETSI code uses a global Overflow flag. In the actual
1634f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber* implementation a pointer to the overflow flag is passed into the function.
1644f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber* This pointer is in turn passed into the basic math functions such as add(),
1654f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber* L_shl(), L_shr(), sub() called by this module.
1664f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1674f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    e = norm_l (L_x) & 0xFFFE;              // get next lower EVEN norm. exp
1684f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    L_x = L_shl (L_x, e);                   // L_x is normalized to [0.25..1)
1694f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    *exp = e;                               // return 2*exponent (or Q1)
1704f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1714f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    L_x = L_shr (L_x, 9);
1724f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    i = extract_h (L_x);                    // Extract b25-b31, 16 <= i <= 63
1734f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber                                                because of normalization
1744f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    L_x = L_shr (L_x, 1);
1754f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    a = extract_l (L_x);                    // Extract b10-b24
1764f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    a = a & (Word16) 0x7fff;
1774f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1784f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    i = sub (i, 16);                        // 0 <= i <= 47
1794f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1804f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    L_y = L_deposit_h (table[i]);           // table[i] << 16
1814f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    tmp = sub (table[i], table[i + 1]);     // table[i] - table[i+1])
1824f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    L_y = L_msu (L_y, tmp, a);              // L_y -= tmp*a*2
1834f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1844f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    return (L_y);
1854f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber}
1864f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1874f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber------------------------------------------------------------------------------
1884f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber RESOURCES USED [optional]
1894f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1904f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber When the code is written for a specific target processor the
1914f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber the resources used should be documented below.
1924f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1934f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber HEAP MEMORY USED: x bytes
1944f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1954f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber STACK MEMORY USED: x bytes
1964f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
1974f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber CLOCK CYCLES: (cycle count equation for this function) + (variable
1984f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber                used to represent cycle count for each subroutine
1994f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber                called)
2004f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber     where: (cycle count variable) = cycle count for [subroutine
2014f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber                                     name]
2024f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2034f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber------------------------------------------------------------------------------
2044f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber CAUTION [optional]
2054f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber [State any special notes, constraints or cautions for users of this function]
2064f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2074f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber------------------------------------------------------------------------------
2084f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber*/
2094f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2104f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas HuberWord32 sqrt_l_exp(      /* o : output value,                          Q31 */
2114f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    Word32 L_x,         /* i : input value,                           Q31 */
2124f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    Word16 *pExp,       /* o : right shift to be applied to result,   Q1  */
2134f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    Flag   *pOverflow   /* i : pointer to overflow flag */
2144f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber)
2154f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2164f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber{
2174f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    Word16 e;
2184f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    Word16 i;
2194f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    Word16 a;
2204f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    Word16 tmp;
2214f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    Word32 L_y;
2224f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2234f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    /*
2244f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber          y = sqrt(x)
2254f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber          x = f * 2^-e,   0.5 <= f < 1   (normalization)
2264f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber          y = sqrt(f) * 2^(-e/2)
2274f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber          a) e = 2k   --> y = sqrt(f)   * 2^-k  (k = e div 2,
2284f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber                                                 0.707 <= sqrt(f) < 1)
2294f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber          b) e = 2k+1 --> y = sqrt(f/2) * 2^-k  (k = e div 2,
2304f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber                                                 0.5 <= sqrt(f/2) < 0.707)
2314f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber     */
2324f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2334f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    if (L_x <= (Word32) 0)
2344f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    {
2354f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber        *pExp = 0;
2364f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber        return (Word32) 0;
2374f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    }
2384f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2394f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    e = norm_l(L_x) & 0xFFFE;               /* get next lower EVEN norm. exp  */
2404f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    L_x = L_shl(L_x, e, pOverflow);         /* L_x is normalized to [0.25..1) */
2414f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    *pExp = e;                              /* return 2*exponent (or Q1)      */
2424f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2434f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    L_x >>= 10;
2444f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    i = (Word16)(L_x >> 15) & 63;            /* Extract b25-b31, 16<= i <=63  */
2454f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    /* because of normalization       */
2464f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2474f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    a = (Word16)(L_x);                      /* Extract b10-b24 */
2484f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    a &= (Word16) 0x7fff;
2494f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2504f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    if (i > 15)
2514f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    {
2524f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber        i -= 16;                              /* 0 <= i <= 47                   */
2534f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    }
2544f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2554f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    L_y = L_deposit_h(sqrt_l_tbl[i]);       /* sqrt_l_tbl[i] << 16            */
2564f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2574f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    /* sqrt_l_tbl[i] - sqrt_l_tbl[i+1]) */
2584f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    tmp = sub(sqrt_l_tbl[i], sqrt_l_tbl[i + 1], pOverflow);
2594f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2604f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    L_y = L_msu(L_y, tmp, a, pOverflow);    /* L_y -= tmp*a*2                 */
2614f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2624f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    /* L_y = L_shr (L_y, *exp); */          /* denormalization done by caller */
2634f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
2644f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber    return (L_y);
2654f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber}
2664f1efc098cb5791c3e9f483f2af84aef70d2d0a0Andreas Huber
267