1 /* ----------------------------------------------------------------------
2 * Project: CMSIS DSP Library
3 * Title: arm_std_q15.c
4 * Description: Standard deviation of an array of Q15 vector
5 *
6 * $Date: 23 April 2021
7 * $Revision: V1.9.0
8 *
9 * Target Processor: Cortex-M and Cortex-A cores
10 * -------------------------------------------------------------------- */
11 /*
12 * Copyright (C) 2010-2021 ARM Limited or its affiliates. All rights reserved.
13 *
14 * SPDX-License-Identifier: Apache-2.0
15 *
16 * Licensed under the Apache License, Version 2.0 (the License); you may
17 * not use this file except in compliance with the License.
18 * You may obtain a copy of the License at
19 *
20 * www.apache.org/licenses/LICENSE-2.0
21 *
22 * Unless required by applicable law or agreed to in writing, software
23 * distributed under the License is distributed on an AS IS BASIS, WITHOUT
24 * WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
25 * See the License for the specific language governing permissions and
26 * limitations under the License.
27 */
28
29 #include "dsp/statistics_functions.h"
30
31 /**
32 @ingroup groupStats
33 */
34
35 /**
36 @addtogroup STD
37 @{
38 */
39
40 /**
41 @brief Standard deviation of the elements of a Q15 vector.
42 @param[in] pSrc points to the input vector
43 @param[in] blockSize number of samples in input vector
44 @param[out] pResult standard deviation value returned here
45
46 @par Scaling and Overflow Behavior
47 The function is implemented using a 64-bit internal accumulator.
48 The input is represented in 1.15 format.
49 Intermediate multiplication yields a 2.30 format, and this
50 result is added without saturation to a 64-bit accumulator in 34.30 format.
51 With 33 guard bits in the accumulator, there is no risk of overflow, and the
52 full precision of the intermediate multiplication is preserved.
53 Finally, the 34.30 result is truncated to 34.15 format by discarding the lower
54 15 bits, and then saturated to yield a result in 1.15 format.
55 */
56 #if defined(ARM_MATH_MVEI) && !defined(ARM_MATH_AUTOVECTORIZE)
arm_std_q15(const q15_t * pSrc,uint32_t blockSize,q15_t * pResult)57 ARM_DSP_ATTRIBUTE void arm_std_q15(
58 const q15_t * pSrc,
59 uint32_t blockSize,
60 q15_t * pResult)
61 {
62 q15_t var=0;
63
64 arm_var_q15(pSrc, blockSize, &var);
65 arm_sqrt_q15(var,pResult);
66 }
67 #else
arm_std_q15(const q15_t * pSrc,uint32_t blockSize,q15_t * pResult)68 ARM_DSP_ATTRIBUTE void arm_std_q15(
69 const q15_t * pSrc,
70 uint32_t blockSize,
71 q15_t * pResult)
72 {
73 uint32_t blkCnt; /* Loop counter */
74 q31_t sum = 0; /* Accumulator */
75 q31_t meanOfSquares, squareOfMean; /* Square of mean and mean of square */
76 q63_t sumOfSquares = 0; /* Sum of squares */
77 q15_t in; /* Temporary variable to store input value */
78
79 #if defined (ARM_MATH_LOOPUNROLL) && defined (ARM_MATH_DSP)
80 q31_t in32; /* Temporary variable to store input value */
81 #endif
82
83 if (blockSize <= 1U)
84 {
85 *pResult = 0;
86 return;
87 }
88
89 #if defined (ARM_MATH_LOOPUNROLL)
90
91 /* Loop unrolling: Compute 4 outputs at a time */
92 blkCnt = blockSize >> 2U;
93
94 while (blkCnt > 0U)
95 {
96 /* C = A[0] * A[0] + A[1] * A[1] + ... + A[blockSize-1] * A[blockSize-1] */
97 /* C = A[0] + A[1] + ... + A[blockSize-1] */
98
99 /* Compute sum of squares and store result in a temporary variable, sumOfSquares. */
100 /* Compute sum and store result in a temporary variable, sum. */
101 #if defined (ARM_MATH_DSP)
102 in32 = read_q15x2_ia (&pSrc);
103 sumOfSquares = __SMLALD(in32, in32, sumOfSquares);
104 sum += ((in32 << 16U) >> 16U);
105 sum += (in32 >> 16U);
106
107 in32 = read_q15x2_ia (&pSrc);
108 sumOfSquares = __SMLALD(in32, in32, sumOfSquares);
109 sum += ((in32 << 16U) >> 16U);
110 sum += (in32 >> 16U);
111 #else
112 in = *pSrc++;
113 sumOfSquares += (in * in);
114 sum += in;
115
116 in = *pSrc++;
117 sumOfSquares += (in * in);
118 sum += in;
119
120 in = *pSrc++;
121 sumOfSquares += (in * in);
122 sum += in;
123
124 in = *pSrc++;
125 sumOfSquares += (in * in);
126 sum += in;
127 #endif /* #if defined (ARM_MATH_DSP) */
128
129 /* Decrement loop counter */
130 blkCnt--;
131 }
132
133 /* Loop unrolling: Compute remaining outputs */
134 blkCnt = blockSize % 0x4U;
135
136 #else
137
138 /* Initialize blkCnt with number of samples */
139 blkCnt = blockSize;
140
141 #endif /* #if defined (ARM_MATH_LOOPUNROLL) */
142
143 while (blkCnt > 0U)
144 {
145 /* C = A[0] * A[0] + A[1] * A[1] + ... + A[blockSize-1] * A[blockSize-1] */
146 /* C = A[0] + A[1] + ... + A[blockSize-1] */
147
148 in = *pSrc++;
149 /* Compute sum of squares and store result in a temporary variable, sumOfSquares. */
150 sumOfSquares += (in * in);
151 /* Compute sum and store result in a temporary variable, sum. */
152 sum += in;
153
154 /* Decrement loop counter */
155 blkCnt--;
156 }
157
158 /* Compute Mean of squares and store result in a temporary variable, meanOfSquares. */
159 meanOfSquares = (q31_t) (sumOfSquares / (q63_t)(blockSize - 1U));
160
161 /* Compute square of mean */
162 squareOfMean = (q31_t) ((q63_t) sum * sum / (q63_t)(blockSize * (blockSize - 1U)));
163
164 /* mean of squares minus the square of mean. */
165 /* Compute standard deviation and store result in destination */
166 arm_sqrt_q15(__SSAT((meanOfSquares - squareOfMean) >> 15U, 16U), pResult);
167 }
168 #endif /* defined(ARM_MATH_MVEI) */
169
170 /**
171 @} end of STD group
172 */
173
